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ABSTRACT 

An existing one-dimensional model of the annual water balance 
is reviewed. Slight improvements are made in the method of calculating 
the bare soil component of evaporation, and in the way surface retention 
is handled. A natural selection hypothesis, which specifies the equilib- 
rium vegetation density for a given, water limited, climate-soil system, 
is verified through comparisons with observed data and is employed in 
the annual water balance of watersheds in Clinton, Ma. , and Santa Paula, 
Ca. , to estimate effective areal average soil properties. Comparison 
of CDF’s of annual basin yield derived using these soil properties with 
observed CDF's provides excellent verification of the soil-selection 
procedure- This method of parameterization of the land surface should 
be useful with present global circulation models, enabling them to 
account for both the non-linearity in the relationship between soil 
moisture flux and soil moisture concentration, and the variability of 
soil properties from place to place over the earth's surface. 
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Chapter 1 


INTRODUCTION 

In order to increase the accuracy of global climate models, 
a more sophisticated representation of the land surface boundary condi- 
tion is required than that which is presently employed (CARP, 1978). The 
interaction, in particular the water flux, between the atmosphere and 
the soil-vegetation system at this boundary is highly non-linear in 
nature, and is not simply defined. Any attempt to satisfactorily 
account for this non-linearity in a model must incorporate txjo effects 
which are not included in current models: 

1. variability of soil properties and soil moisture 
dynamics from place to place over the earth's surface, and 

2. non-linearity in the relationship between soil moisture 
flux and soil moisture concentration. 

In this work, it is intended to make use of a one-dimensional 
water balance at the land-air interface in order to parameterize the 
climate-soil-vegetation relationship in such a way as to reflect the 
non-linearity and areal variability. 


14 



Chapter 2 


OBJECTIVES 

The specific objectives of this work are twofold: 

The first objective is verification of a vegetal equilibrium 
hypothesis developed by Eagleson (1978f). This hypothesis proposes 
that the natural vegetation density in a watershed will seek, through 
natural selection, an optimal "climax" value at which available soil 
moisture is a maximum. Comparison of a theoretical curve of evapo- 
transpiration versus canopy density based on this hypothesis x^ith 
observed data will provide the necessary check on the accuracy of the 
hypothesis. 

The second objective is establishment of an algorithm for 
estimating effective areal soil properties from observations of 
vegetation density by using the natural selection hypothesis in a 
one-dimensional water balance model. By defining the level of evapo- 
transpiration from soil moisture through observations of the canopy 
cover density, it may be possible, knowing the climate, to determine, 
the soil properties that enable the soil-vegetation system to respond 
at the indicated level. The estimated values of these parameters can 
then be used in the water balance equation to evaluate desired compo- 
nents of the water flux. Verification of the desired algorithm will 
be sought through comparison of computed and observed statistics 
of annual yield. 
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Chapter 3 


REVIEW OF LITERATURE 

Past efforts to model the coupling among physical processes 
of the atmosphere, soil, and vegetation across the land-surface inter- 
face have been largely of two types: 

1. Numerical studies which employ detailed formulations 
of the processes involved. Examples of such studies are those of 
Philip (1957), Sasamorl (1970), Deardorff (1977), and Philip and de Vries 
(1957) . Although these models simulate the system response to climatic 
inputs very well, they usually do so in terms of a large number of 
climate, soil, and vegetal parameters. Due to their complexity and 

the detailed data requirements for their validation,, these studies 
have found little application in general circulation models. 

2. Empirical studies which utilize validated interrelation- 
ships among the principle variables. Because of the ease of their 
application, and negligible programming and data requirements, most 
global climate models use this type of parameterization of the land- 
surface boundary with regards to actual evapo transpiration, average 
soil moisture content, and runoff. 

The primary GCM’s today utilize the approach first introduced 
by Manabe (1969) to parameterize the land surface boundary condition. 

In this approach, the above mentioned parameters are handled in the 
following way. 
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A. Evapo transpiration 

Actual evapotranspiration is related to potential evapo- 
transpiration linearly through the soil moisture and a single soil 
parameter following the work of Budyko (1956) . This parameterization 
is 


e 


T 




s/k, s ^ k 
1 , k < s < 1 


in which 


(3.1) 


e^ = actual rate of evapotranspiration 
Bp - potential rate of evapotranspiration 
s = effective soil moisture concentration 
k = empirical coefficient, 0<k ^ 1 generally assumed 
to be constant everjn^here 


As mentioned above, the only soil parameter appearing in this 
model is the empirical coefficient, k. This representation grossly 
distorts the sensitivity of to s and makes no allowance for the spatial 
variance of this sensitivity due both to soil type and to the presence 
of vegetation. 

More recently, Lettau (1969) and Lettau and Baradas (1973), 
in their "evaporation climatonomy" formulation, refine the water 
balance evapotranspiration term through use of an energy balance. This 
approach seeks theoretical solutions in the form of "response functions" 
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(i.e. evapotranspiration cycles, temperature, etc.) as a physical con- 
sequence of a mathematically defined "forcing function" of the environ- 
mental system. However, parameterlsation is achieved without any 
explicit consideration of the soil and vegetal properties which x<rill 
control the evaporation under all but the most humid conditions. The 
only input parameters linked to the land surface are evaporivity, e*, 
which is a non-dimensional measure of the capacity of land surfaces to 
utilize solar energy for the evaporation of rainfall received in a 
specified time interval, and t* which denotes a characteristic soil 
moisture residence time. Values for these parameters are either 
assumed on the basis of empirical data, or are estimated from a systematic 
classification of watersheds according to morphology, soil structure and 
permeability, vegetation cover, etc. The lumping of all these para- 
meters into a single term in no way fully represents the complex 
interrelationships between the various processes involved in the 
water balance. 

Other studies concerning the evapotranspiration term are 
those of Czarnowski (1964) and Ritchie (1972) , and Ritchie and Burnett 
(1971) . Czarnowski assumes that total evapotranspiration is a sum of 
plant transpiration, evaporation from surface retention, and evaporation 
from soil, and that these values are functions of vegetation density, 
and consequently of climatic factors. He treats the development of 
plant cover as a function of the form 

_ JL 

V 

M = 1 - e “ (3.2) 
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where 


V — precipitation, mm 

= sum of mean daily deficits of air humidity, mm Hg 
Finally, he concludes that evapotranspiration can be expressed as 


V 


V 

m 


M 


r 

|l.l7M + 



or 


V 4 

= J = M 1.17M + 

m /W 


(3.3) 


(3.4) 


where the constants, 1.17 and .4 are determined by a least squares fit 
to empirical data obtained primarily from cultivated agricultural 
lands . 

Ritchie (1972) and Ritchie and Burnett (1971) develop a set 
of empirical functions relating leaf area index and fractional net 
radiation at a soil surface for a row crop to plant evaporation 
efficiency. These equations may be ivuritten 


R 


ns 


-.398L 


Ai 


e 


(3.5) 


— = .70L - .21 (3-6) 

Op Ax 


where 


L»* • 

Ax 


leaf area index 
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net radiation reaching soil surface 


= net radiation above plant canopy 

Relating leaf area index to canopy density using Equation (3.5) and 
the assumption that 


R 

ns 


1 - M 


(3.7) 


gives 


J =» 



.70 


In(l-M) 

-.398 


1/2 


.21 


(3.8) 


Again, the constants appearing in the above equations are determined 
from the method of least squares. 

, l-Jhile both of the above formulations are attempting to relate 
evapotranspiration to more physically significant parameters, there is 
little inclusion of the actual physics. Since a linear regression is 
performed to obtain the above equations, there is a lack of generality 
and understanding of the sensitivity to other parameters besides vege- 
tation 


B. Soil Moisture 

The change in the average soil moisture concentration is 
determined from a water balance relation written 
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" - "t - ■"s - 


<3.9) 


nh 


9s 

9t 


In which 

h = thickness of surface layer 

i = intensity of rainfall 

Yg = intensity of surface runoff 

Yg = intensity of percolation of water out of the surface 
layer to groundwater 

The product, nh, represents the maximura water content of 
the surface layer and is assigned a value which is conmion for all 
soil surfaces. 

C. Runoff 

Runoff, as written in Equation (3.9; consists of two different 
components. Surface runoff is regulated hy the infiltration of rainfall 
and additions to soil moisture. Groundwater runoff is governed by the 
state of soil moisture concentration. All global climate models use 
highly simplified empirical formulae which Itimp these two dynamically 
different runoff-generating processes into total yield relations of 
the form 


Y = Y(i,e^, s) 
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These relations include one or more coefficients which may incorporate 
spatial variability, but there is no physical basis for their selection 
without natural yield measurements. 

The models referred to in the preceding paragraphs include 
those of Arkawa (1972, U.C.L.A.); Somerville et al. (1979, G.I.S.S.); 
Gates and Schlesinger (1977, Rand-‘0. S .U. ) ; Sellers (1973, Arizona); 
and Corby et al. (1978, B.M.O.). In all of these models, there is 
no use of the present high level of physical understanding of the 
natural processes involved to develop a generalized, accurate repre- 
sentation of the land-surface interface. 

Eagleson (1978a, b, c, d, e,f,g) , has developed a generalized 
water balance based upon simplified physics of the component processes. 
The development is sufficiently rigorous to capture the essential 
system dynamics yet simple enough to permit analytical solution. 

The model produces valuable insights into the interactive role of 
soil moisture in the determination of climate. Foremost in this 
development is the accounting for the areal variability of soil pro- 
perties over the earth's surface and the reflection of the inherent 
non-linearity in the relationships betireen soil moisture concentration 
and the interfacial moisture fluxes . This model will be presented 
and utilized in the following chapters to attain the objectives 
Stated in Chapter 2. 
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Chapter 4 


THEORETICAL BACKGROUND 

4.1 The Water Balance 

The major source of theoretical background used here is the 
work of Eagleson (1978a, b, c, d, e, f ,g) . In these papers, a one- 
dimensionaL water balance based on soil moisture dynamics and 
statistics of climatic data is derived. This water balance, expressed 
in terms of annual expected values, may be represented as 


E[Ep ] = E[E^ ] + E[Rg ] + E[Rg ] 

A A A A 


and 


E[Y^] = E[Rg ] + E[Rg ] 


where 


E[ ] 



expected value of [ ] 

annual precipitation 
annual evapotranspiration 
annual surface runoff 

annual groundwater runoff 

annual yield 


(4.1) 


(4.2) 


23 



An analytic expression is obtained for Equation (4.1) by deriving the 
individual components through the use of derived probability distribu- 
tions and one-dimensional dynamic equations approximating the physics 
of the separate soil moisture fluxes. These expressions are then 
introduced into Equation (4.1) to produce the equation for the (soil 
moisture) water balance 

HL r(a + 1)0“^) = 

A 

E[E^ ] J(E,M,k ,h ) - E[E ] + m K(l) s Tt<r 
P. ^ ’ V* o r. T o 

A A 


for 


E[E ] / EL < e r(0 + l)a ^ (4.3a) 

A A 

(The term to the left of the equal sign is infiltration, the first 
term to the right, total evapotranspiration, the second is evaporation 
from surface retention, and the last two terms are groundwater runoff 
[the first is groundwater recharge and the last is groundwater loss]). 
Otherrvfise, 


mp - E[Ep ] J(E,M,k^,h^) + m^K(l)s^‘' - Tw (4.3b) 

A A 

In the above, 

E = annual surface retention 
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E = average annual potential evapo transpiration 
J = evapotranspiration efficiency 
G = gravitational infiltration parameter 
a ~ capillary infiltration parameter 
E = evapotranspiration parameter 
M = vegetation canopy density 
k.^ = plant transpiration coefficient 
= mean length of rainy season 
- surface retention capacity 
= average annual soil moisture 
K(l) = saturated hydraulic conductivity 
T = 1 year, seconds 

w - apparent velocity of capillary rise 




mean annual precipitation 


pore disconnectedness index 


It will he helpful and important to review the development of 
the expressions for evapotranspiration and surface runoff, and to 
present an alternative approach for the former, and a slightly different 
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interpretation for the latter. 

4,2 Evapotranspiration 

The expected value of annual evapotranspiration is derived 
(Eagleson, 1978d) by calculating bare soil evaporation and vegetal 
transpiration for an interstorm period as functions of properties of 
the storm sequence, the surface, the soil, and the average rate of 
potential evapotranspiration, using observed distributions of the 
random climatic variables » and averaging over the rainy season. The 
bare soil evaporation and plant transpiration are determined by con- 
sidering the vertical flux of moisture in a soil column. In Figure 
4.1, the modeled column of soil and the different moisture fluxes 
are sketched. In this figure 

f = hare soil exfiltration rate 
e 

M = vegetation canopy density 

e^ = vegetation transpiration rate 

K(0^) = effective hydraulic conductivity at long-term average 
soil moisture 

It is assumed here that 

1. Soil moisture throughout the surface boundary layer is 
spatially uniform at the start of each interstorra period 
at the long-term average value, s = 

2. The medium is semi- infinite; 


26 



fe(t) 

t t t 

— I-M 



V 11 V I 


K (0J 


Figure 4. 


SCHEtlATIC REPRESENTATION OF 
DURING AN INTERS'] 


2 



3. The vegetation is distributed uniformly, and its roots 
extend uniformly into the entire volume of the soil in 
the surface boundary layer. This implicitly assumes that 
the plant species have adapted by natural selection to 

a density and root structure which is in balance i-riLth 
the available soil moisture; 

4. The rate of moisture extraction by the roots is in 
equilibrium with the transpiration rate by the leaves, 
thus forming a uniformly distributed sink for soil 
moisture of strength, Me^- 


Following the work of Philip (1969) , Eagleson writes the 
total decrease in soil moisture during infiltration: 


(0^ - 0)dz = 


zd9 = F^(t) + [K(0^) + Me^I t 


(4.4) 


where F (t) is the cumulative exfiltration in centimeters. 

The integral on the left-hand side is evaluated in the manner 
of Philip (1960) . Assuming a vertical flow passage of constant 
cross-section, the exfiltration rate is found to be 

_1 

fe(t) =1 S^t ^ - |lK(0^) - K(0^)] - Me^ (4.5) 


Note that this neglects the restriction, by vegetation canopy density, 
of the bare soil area through which exfiltration occurs. Further 
simplification and analysis result in the exfiltration capacity: 
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(4.B) 


1 2 

•r S t -'Me 

2 e V 


where is the exf iltratioa sorptivity. 

A typical interstorm period, and the relationship between the 
various rates of evaporation and time for bare soil is illustrated 
in Figure 4.2. In this figure, is the potential rate of bare 

soil evaporation, which is considered a constant. The times, t^ and 
tg, are evaluated by assuming 


e = f 
p e 


(4.7a) 


when 


f*(Ve) 



(e t ) 

^ P o 


(4.7b) 


and that 


f^ = 0 (4.8a) 

when 

t = t^ (4.8b) 

respectively. Exfiltration capacity, and the times t^ and t^, are 
then used by Eagleson along with the relationships represented in 
Figure 4.2 to determine total evapo transpiration, E^. To do this. 
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Figure 4.2 
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E.J from a unit land surface is proportioned according to 


^ (l-M)E + tlE 
T s V 


(4.9) 


in which 


Eg = bare soil evaporation from soil moisture plus eva- 
poration from soil surface retention 
= evapo transpiration from vegetation plus evaporation 
from plant surface retention. 


It is not necessary to present the development of E[Eg] here, 
which is done by calculating the volume under the solid line in Figure 
4.2, multiplying by the joint probability distribution of storm depth 
and time between storms, and integrating over the regions sho\m in 
Figure 4.3. What is important to note is the previously mentioned 
approximation made in the development of the bare soil exfiltration 


capacity. The expression obtained for E[E ] , bare soil evapotrans- 

S • 

3 

piration fof one inter storm period, from the above procedure is 



I'-ytKjXh^] 


-BE-6h Is 

° P T /2 
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+ [Hk^ + (2C)^^^ E - w/e^] - [Mk^ 4- (2B)^^^E ~ w/e^] 


Here 


(4.10) 


and 


B = 


1-M 


1+Mk -w/e 

V p 


M k +( 1-M) w/e 

V p 

2(1+Mk -w/e y 
V p 


(4.11) 


C = (Mk - w/e ) ^ 
2 V p 


(4.12) 


Also, 


B = 
X = 

K = 


h = 
o 


reciprocal of mean time between storms 
parameter of Gamma distribution of storm depth 
parameter of Gamma distribution of storm depth 
surface retention capacity 
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Upon studying Figure 4.1 and Equations (4.4) and (4.6) 

it can be seen that the term, f , is defined as the exfiltration rate 

e 

for bare soil, and as the total volume of moisture exfiltrated from 

the soil column across the bare soil surface. The rate, f , is 

e 

obtained by differentiating the volume, F^, x-^ith respect to time. 

The result of the differentiation leaves f^ multiplied by the area 

of bare soil. Thus, in the two-dimensional problem which includes 

the presence of vegetation, f should be multiplied by the terra, 1-M, 

e 

to account for the fact that only a fraction of the land surface, the 
bare soil fraction, is exfiltrating at this rate. Equation (4.6) 
should be rewritten as 

(l-M)f = ^ S t ^ - Me ■ (4.13) 

e2e V ' 


The new form for the expected value of bare soil evaporation, E^, may 
then be evaluated in terms of this altered expression for exfiltration. 


The new expression for E[E ] is 

3 


= f ^ 
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(4.14) 


where 

B = 1/(1+Mk^) (4.15) 

C + 1/2 (Mk -w/e ) ^ (4.16) 

V p 

is obtained in the same way as before; by multiplying the bare soil 
term by (1-M), and the vegetation term by M. The result of this altera- 
tion on the expected value of annual basin evapotranspiration will be 
presented in a later section. Although this approach may seem more 
accurate than the original, its use will create other, and possibly 
greater problems . Attempting to expand the problem into two dimensions 
at this point will cause some inconsistencies concerning evaluation of 

the Philip exfiltration equation. This equation is 

6 

o 

zd0 = A + A^t + + ... (4.17) 

Since this was developed for a one-dimensional formulation, the 
expressions obtained for the constants, A^, on the right hand side will 
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not necessarily apply to- the two-dimensional situation. This can be 
seen by noting that in Equation (4.6), as M approaches the value of 1, 
the right hand side does not go to zero, as it should for a fully 
vegetated surface x^here there x-rould be no bare soil exfiltration of 
soil moisture. So, although there are certain misgivings about 
Eagleson’s original derivation, the alternate approach presented above 
may involve more serious inaccuracies. Hoxjever, for areas with a large 
vegetal canopy density, xAere the effect of the vegetation on bare soil 
exfiltration is large, this approach may come closer to reality than 
the previous one. 

4.3 Surface Runoff 

To derive the probability of storm surface runoff, Eagleson 
(1978) integrates the difference between rainfall intensity and the 
Philip infiltration equation over the duration of a rainstorm. Infil- 
tration is assumed to occur uniformly over both bare soil and vegetated 
portions of the surface. Illustrated in Figure 4.4 is a sequence of 
surface states beginning from t = 0 at the start of the rainfall period. 
In this figure 

h = surface retention capacity 

0 

1 = rainfall Intensity 

f^ = infiltration capacity 

t = storm duration 
r 

= gravitational infiltration rate as modified by capillary 
rise from the xjater table 


36 




SURFACE RUNOFF GENERATION DURING TYPICAL 
STORM (tr>to) 


Figure 4.4 



Initially, there is a withdrawal of rainfall to satisfy surface retention. 

If t^ is greater than h^/i, as shown, this surface retention reaches its 

capacity, h^. If t^ j< h^/i, there would be no infiltration or runoff, 

and the surface retention would equal the storm depth, h. For the case 

illustrated in Figure 4.4, however, infiltration will begin at time 

t = h /i. From this time until t = h /i + t , when f , = i, infiltration 
o o o X 

will take place at the rainfall rate, i. After this time, the capacity 

of the soil to infiltrate moisture is no longer larger than the rainfall 

intensity, and the excess is represented by the shaded area of the 

figure. Rainfall excess ^ R , is then generated until time t = t . The 

s * r 

3 

expected value of the rainfall excess is obtained in a manner similar to 
that of the evapo transpiration. 

A question may be raised relating to the handling of the 
surface retention. In his development, Eagleson argues that the surface 
retention must be subtracted from the rainfall excess, since it is 
moisture that is not infiltrated into the soil. The expression he 
obtained for the expected value of annual rainfall excess is then 

E[R* 3 = in r(a + l)a“^3 (4.18) 

®A A 

in which 

R. = annual rainfall excess 
^A 

The expected value of annual surface retention, E , is then 
subtracted from this to get the annual surface runoff. This charges the 
entire annual surface retention against those events producing rainfall 
excess, however. 
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CLAY 

CLAY-LOAM 

SILTY-LOAJI 

SANDY-LOAM 

k(l) 

1.0x10”^° 

-10 

2.8x10 

1.2x10"^ 

2.5x10"^ 

n 

.45 

.35 

.35 

.25 

c 

12 

10 

6 

4 


Table 4.1 

REPRESENTATIVE SOIL PROPERTIES 
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A slightly different interpretation of this rainfall excess 
results from a closer examination of Figure 4.4 and the shaded area 
therein. It is kno\m, and mentioned above, that surface retention 
must be satisfied before any infiltration can occur. From this, it 
seems necessary to subtract the surface retention from the beginning 
of the rainfall period, as indicated in the figure, rather than from 
the rainfall excess at the end. The volume represented by the shaded 
area would then be equal to the surface runoff, and not surface runoff 
plus surface retention. The resulting water balance equation then 
becomes 

% (1 - = E[E ] J(E,M,k h ) + m KCDs*^ " Tw (4.19) 

V o T o 

This alternative procedure will increase the calculated value of sur~ 
face runoff, and decrease the amount of moisture calculated as infiltra- 
tion. The effect of this difference on the CDF of annual yield will 
be discussed in Chapter 6. 

4.4 Vegetal Equilibrium Hypothesis 

From examining the role of vegetation canopy density in the 

average annual water balance, Eagleson (1978f) observed that for a given 

set of climate and soil parameters and for a given k^. Equation (4.3) 

defines as a function of M. This relationship is illustrated in 

Figure 4.5 by using four sets of representative soil properties, listed 

in Table 4.1, and the conditions ~ for the climates 

A 

of Clinton, Mass, and Santa Paula, Calif. It can be seen that there 
exists a particular value of M = M for each climate-soil combination at 
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Figure 4.5 



FOR TYPICAL SOILS (k^ = 1 AND u/e^ « 1) 


41 





which s is a maximum. This point of maximum s corresponds to 
o o 

maximum surface and groundwater runoff, which means, for fixed pre- 
cipitation, that there is minimum evapotranspiration from soil moisture. 
Thus, at M = M^, it is expected that the term representing evapo- 
transpiration from soil moisture 


E[E ]J(E,M,k^,ho) - E[E ] 

Pa a 


(4.20) 


will be a minimum for a given climate-soil combination. This minimiza- 
tion is seen in Figure 4.6 for the same information as that used in 
Figure 4.5. Note that in Santa Paula, the clay and clay loam soils 
cannot absorb enough water to produce canopy densities greater than 0.4 
and 0.8, respectively, as long as = 1. 

The numerical value of is a matter of some controversy. 
Linacre, et al. (1970), report values of k^ for water plants which range 
from .6 to 2.5 depending upon species. Slatyer (1967, p. 53) states 
that the value of k^ can be greater than one since total evapotranspira- 
tion from a plant community, per unit land area, may exceed that from a 
similar area of bare wet soil due to the larger actual evaporating 
surface area. Kramer (1969, p. 338) however, states that evaporation from 
a plant community never exceeds that from a similar area of wet soil. Por 
the present, k^ will not be allowed to exceed one. 

•From observations of the relationships presented above, 

Eagleson (1978f) develops the vegetal equilibrium hypothesis mentioned 
in Chapter 2. In the light of the above arguments, this hypothesis says 
that natural vegetal systems of given species will develop a canopy density 
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which produces minimum stress under local climatic conditions. A 
necessary condition for minimum stress is that soil moisture take on. 
the largest possible value. Thus, by using this hypothesis, the given 
climate, soil, and plant coefficient determine the equilibrium canopy 
density, M == M^, through the water balance equation, where the soil 
moisture is maximum or, equivalently, where the soil moisture evapo- 
transpiration is a minimum. Figure 4.7 illustrates the relationship 
between the dimensionless evapotranspiration parameter, E, and the 
dimensionless evapotranspiration function, J(E, kv) , for the equilibrium 
condition, M = This plot is obtained by minimizing evapotranspiration 

from soil moisture for a given k^, (k^ = 1 in this case) , and E using 
Eagleson's constant soil column cross-section assumption. The expression 
for E is 


E = 


2Bn KCD "FCD 
_____ 

TTme 

P 


e d+2 
®o 


(4.21) 


in which 

3 = reciprocal of mean time between storms 
n = porosity 

t!*(l) = saturated soil matrix potential 
m = soil pore size distribution index 
d = soil diffusivity index 
(j)g = dimensionless desorption diffusivity 

The other terms have been previously defined. 

Also shoxm in Figure 4.7 is the vs.E relationship for the 

equilibrium condition, P. = m^ . As Pa varies from mp , E and thus, 

^ ^A * 
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J(E, k^) will change accordingly, while to the first approximation, 

M will remain constant at M . 

o 

Eagleson (1978f) performs an asymptotic behavior analysis of 
the evapotranspiration function to gain insight into the meaning of the 
parameter, E. The evapotranspiration asjnap totes sho^-m in Figure 4.7 
are thereby determined. The intersection of these ti?o asymptotes occurs 
at E = 2 /tt, which separates soil controlled from climate-controlled 
evapotranspiration (Eagleson, 1978d) . Thus, low values of E correspond 
to relatively dry, warm climates, while larger values indicate humid 
climates. As can be seen from Figure 4.7, low values of occur for 
low E values, and vegetation densities approaching 1 correspond to a 
large E. 

It can now be seen that observations of canopy cover will 

provide a key to determining the effective properties of a soil for a 

given climate. By using the vegetal equilibrium hypothesis in reverse, 

observations of M may be used in the i/ater balance to obtain information 
o 

about the soil if the climatic variables are known. 

Figure 4.8, whicli is a plot of J vs. can be obtained 
directly from the information in Figure 4.7. Thus, from observations of 
vegetation density, the evapotranspiration efficiency, J, can be deter- 
mined. To assure the generality of this relationship, the sensitivity of 
J to its independent parameters is studied. From the expression obtained 
by Eagleson (1978d), the primary parameters other than E and M are: 

k = plant coefficient 
V 

Xh accounts for storms which do not fill retention capacity 
o 

3h /e measures effect of surface retention on exfiltration 
o p 

K = parameter of Gamma distribution of storm depth 






Figure 4.8 

TON FMCTION FOR NATURAL SYSTEiiS 

= 0.0, K = . 5 , w/ « 1) 
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The plots obtained by varying these parameters over their reported 
ranges are presented in Figures 4.9, 4.10 and 4.11. Figure 4.9, which 
holds fc and constant at .5 and 1, respectively, illustrates the 
insensitivity of the eyapotranspiration function to changes in Ih^ as 
compared to 3h /e . By holding Bh /e equal to . 1 and k equal to 

op op 'V 

1, K is varied in Figure 4.10. As can be seen, the function changes 

infinitesimally with changes in k. In Figure 4.11, the tx<ro variables 

K and Xh are held constant at median values, and 3h /e and k are 
0 o p V 

allowed to vary. From this analysis, the evapotranspiration function 

is sho\ra. to be most sensitive to the two parameters , Bh /e and k . 

op V 

Also shown in Figure 4.11 as dashed lines are the curves obtained 

using the alternate formulation of evapotranspiration. Equation (4-14), 

developed in Section 4.2. In review, this expression was developed 

by accounting for the effect of the vegetated fraction of the soil 

column surface on the vertical flux of the exfiltrating soil moisture 

in Equation (4.6). Expanding the Philip exfiltration equation, which 

was developed for the one-dimensional case of a constant cross-section, 

to two dimensions introduced an inconsistency with the results Philip 

* 

obtained as explained in Section 4.2. By multiplying the terra f^ in 
Equation (4.6) by (1-M), and not adjusting the terms on the right hand 
side of the expression, an infinite exfiltration capacity is obtained 
for the case when M = 1. Although the terra (1-H) appears in the 
denominator of several components of the equation (4 . 14) for bare soil 
storm exfiltration volume, an infinite result is not obtained since the 
total volume of bare soil exfiltration, E^, is weighted by (1-M) in 
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Figure 4.9 


SENSITIVITY OF EVAPOTRANSPIRATION FUNCTION TO 

(k = 1, K = .5, w/e « 1) 

V P 





AND Ah 

o 
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Figure 4.10 

SENSITIVITY OF EVAPOTRANSPIRATION FUNCTION TO ic 

(Bh /e - .1, k = 1, w/e « 1) 
op V p 
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Figure 4.11 

SENSITIVITY OF EVAPOTEMSPIRATION FUNCTION OF 3h /e = ,1 AND k 

_ op V 

(k = ,5, AN = 0.5, w/e « 1) 

’ o P 
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Equation (4.9). The net result of this modification is' to raise the 
total evapotranspiration for a given value of M since the rate has 
been increased. This is seen in Figure 4.11 where the dashed lines 
are plotted above the corresponding solid lines for the same vegeta- 
tion density. The main problem with this approach, as mentioned 
above, is that the terms on the right hand side of Equation (4-6) 
do not identically go to zero as M approaches 1. If the necessary 
corrections were known, the result would be a reduction in the bare 
soil exfiltration capacity for each value of M. This would lower 
the dashed lines of Figure 4.11. The actual function may therefore 
lie somewhere between these two sets of curves. With this in mind, 
these plots will be used in the following chapters to study the 
validity of the vegetal equilibrium hypothesis, and to determine its 
utility in estimating the effective average areal soil properties of 
a natural watershed. 
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Chapter 5 


METHOD OF APPROACH 


5.1 Vegetal Equilibrium Hypothesis 

Verification of the vegetal equllihriiim hypothesis presented 
in Chapter 4 is the first objective of this work. This can be accomplished 
through comparisons of actual data (from watersheds representing various 
types of climates) with the hypothesized relationship illustrated in 
Figure 4.11. In review, the hypothesis states that the vegetation 
denisty seeks that value, M^, which maximizes soil moisture. This value 
maximizes water yield and thus, for a given climate and soil, minimizes 
evapotranspiration from soil moisture. Minimum evapotranspiration 
can be translated into a value of evaporation efficiency, J. (i.e. the 
ratio of actual to potential evapotranspiration) leading to the 
relationships previously presented in Figure 4.11. 

The average annual water balance is presented by Eagleson 

(1978a) as 

E[P^]-E[EtJ=E[RsJ+E[RgJ=E[Y^] ( 5 ^^^ 

which states that average annual precipitation minus average annual 
evapotranspiration will equal the average annual basin yield which is 
composed of surface runoff plus groundwater runoff. VTien analyzing 
catchment data to calculate average annual actual evapotranspiration, 
mean annual basin yield is subtracted from mean annual precipitation. 
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Although relatively accurate annual precipitation data are readily 
obtained from station records, yield information is only available 
in the form of streamflow records. Therefore, it is necessary that 
the total yield of the catchment studied appear as streamflow. This 
means that the entire groundwater component of yield must be influent 
to the stream channel upstream of the basin mouth. Under such condi- 
tions, most closely approached in humid climates, the total evapo- 
transpiration is equal to precipitation minus streamflow. This 

restriction may lead to overestimating actual evapotranspiration if 

> 

there are losses of yield to ungaged groundwater, or underestimation 
if there is contribution to streamflow of groundwater from adjacent 
watersheds . 

Potential evapotranspiration is estimated by using the 
modified Penman equation (Penman, 1948) . The form used here is the 
combination form as presented by Eagleson (1977) 

1 = 

P Pg (1 +yM) ^ ^ ^ 

in which 

= average rate of insolation 

< 5 ^ = average rate of net outgoing long wave radiation 
H = average residual sensible heat flux 
A = shortwave albedo of surface 

^ mass density of evaporating water 
=• latent heat of vaporization 

y/A = atmospheric parameter, a function of atmospheric temperature 
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The above parameters are calculated or estimated as follows: 

= q^Cfj)); from Figure 5.1, where latitude 

A = A (surface structure); from Table 5.1 

1/(1 + Y/A) = f(T^); from Figure 5.2, where T^ = average 
annual temperature 

L = 597 cal/gr 

Pg = 1 gt/cm^ 

q^ = (1 - .8N)[.245 - .45 X 10"^°T^1 
H = q^(.25 +1/(1 - S)) 

The necessary climatic variables are available from U.S. 
Weather Bureau publications. They must be averaged over the rainy 
season which is assumed to be identical with the vegetation growing 
season. 

Equation (5.2) gives the average potential evapo transpiration 
rate. The total potential volume is obtained by multiplying e^ by the 
season length as determined from monthly rainfall records. 

With actual and potential evapotranspiration known, the only 
remaining variables needed for comparison with the hypothesis are the 
vegetation species (to obtain k^) and the canopy density. The canopy 
density is estimated either from aerial photographs, from personal 
observation, or from literature available for the catchment studied. 

In this work, no photographs were available, and it was possible to 
estimate only ranges of density from the information in the literature, 
depending upon each author's method of measurement and interpretation. 

As a result of this uncertainty regarding the type and canopy density of 
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Figure 5.1 



INSOLATION AT EARTH 

‘ S 

SURFACE 

(from Ref. (15)) 






Table 5-1 


Albedo of Natural Surfaces 
(from Ref. [17]) 


Surface 

Albedo, A 

Surface 

Albedo, A 

Water 

0.03-0.40 

Spring wheat 

0.10-0.25 

Black, dry soil 

0.14 

Winter wheat 

■ 0.16-0.23 

Black, moist soil 

0.08 

Winter rye 

0.18-0.23 

Gray, dry soil 

0.25-0.30 

High, dense grass 

0.18-0.20 

Gray, moist soil 

0.10-0.12 

Green grass 

0.26 

Blue, dry loam 

0.23 

Grass dried in sun 

0.19 

Blue, moist loam 

0.16 

Tops of oak 

0.18 

Desert loam 

0.29-0.31 

Tops ^of pine 

0.14 

Yellow sand 

0.35 

Tops of fir 

0.10 

White sand 

0.34-0.40 

Cotton 

0.20-0.22 

River sand 

0.43 

Rice field 

0.12 

Bright, fine sand 

0.37 

Lettuce 

0.22 

Rock 

0.12-0.15 

Beets 

0.18 

Densely urbanized 


Potatoes 

0.19 

areas 

0.15-0.25 

Heather 

0.10 

Snow 

0.40-0.85 



Sea ice 

0.36-0.50 
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Figure 5.2 

LONG-TERM AVERAGE POTENTIAL EVAPOTRANSPIRATION 

(from Ref, (15)) 




the vegetation, the other variable which is a function of the surface 
structure, albedo, is subject to error as well. Therefore, the 
catchments studied will be plotted on Figure 4.11 in terms of the 
expected range of both J and 

5.2 Estimation of Effective Soil Properties 

The second goal of this work is to use the hypothesized 
relationship between vegetation canopy density and evapotranspiration 
to estimate effective average areal properties of the soil. 

Three types of parameters are, considered: climate, soil 

and vegetation. The climatic and vegetal properties are easily obtained 
from observations; this leaves the four soil parameters, s^, k(l) , n, 
and c to be determined from the derived relationships between climate, 
soil and vegetation. 

The range of values of the porosity, n, is kno\na to be quite 
small, from .25 to about .45, and does not have a large effect on solu- 
tions of the water .balance equation. Assuming a value for n leaves the 
soil moisture, intrinsic permeability, and pore disconnectedness index 
as unknowns. To solve for these variables, three equations or relation- 
ships are needed which incorporate the vegetation and climate as well. 
The first relationship is the water balance. Equation (4.3), which 
expresses the soil moisture, s^, as an implicit function of the climate, 
vegetation and soil. The vegetal equilibrium hypothesis provides the 
second relationship between the same three parameters. The third 
expression used is a rather weakly correlated regression between the 
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SATURATED PERI-IEABILITY VS. PORE SIZE DISTRIBUTION INDEX 
(from Eagleson, Personal Communication) 





intrinsic permeability of the soil, k(l), and its pore size distri- 
bution index, m, presented in Figure 5.3. This expression is 
(Eagleson, Personal Communication) 

'5.3) 

where 

m - 2/(c-3) (5.4) 

The coefficient of determination of this regression is small 
due to the extreme variability of these parameters in nature. The effect ‘ 
of this regression equation on the derived CDF of annual yield will be 
observed in Chapter 6. 

In order to explain the procedure followed in the estimation 
of soil properties, it is necessary to present mathematically the water 
balance and the vegetal equilibrium hypothesis. The mean annual water 
balance. Equation (4.3), is again 


jn (1 _ e ^ ^‘^r(a+l)cy ^) = E[E ] J(E,M,k ,h ) - E[E ] +m^ K(l)s^ - 


for 


E[E ]/nL < r(o + l)a”^ 

-^A 


Otherwise, 


(5.5a) 


nn = E[E j J(E, M, k , h ) + m K(l) s*" - Tw (5.5b) 

v’ o' T ^ ^ o 

If the interpretation of surface runoff developed in Section 
4.3 is used, the above equation becomes 
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lOp (1 - 

A 


e“G~2<J + l)a“^) = E[E ] J(E,M,k ,h ) + m K(l)s^ - T^^.‘ ' 

Va VO T O 

(5.6) 

Although the components and symbols have been defined earlier. 


their full expressions have not all been stated. In the above equations 


G = aK(l)[(l + s )^/2 - w/K(l)] 


(5.7) 


■e 


nTi^ K(l) 'F(l)(l - s^)^ 6.(d, 


(1 “ s^)^- <^^(d, sW 

TT 6 m 


(5.8) 


E[E ] = m m [1 - M(1 - k )] i 

Pa ^ V T 


(5.9) 


K(l) = k(l) y /p 
'w w 


(5.10) 


where 


"^w k(l) <{)(m)-^ 


in which 


a = reciprocal of mean storm intensity 
q = reciprocal of mean storm depth 
d = 2 + 1/m 


6 = reciprocal of mean storm duration 


m = mean number of storms 
V 


m = mean time between storms 
b 


(5.11) 


= specific weight of water 

a = surface tension of water 
w 

= viscosity ^of water 


. , y, , -„0,66 + 0.55/m + 0.14/m 

9 Cm) = pore shape parameter = lU 


62 



The vegetal equilibritim hypothesis states that 


] 

A 

9M 


= 0 


at II = M 

o 


where, as mentioned before. 


^ “ evapotranspiration from soil moisture; 
A 


which is 


E[E ] - J(E, M, k h ) E[E ] - E[E^ ] 
A ^A A 


with 


> - f • ’“v 


(1 


- M)|l - 


-3h /e r[K, Xh ] 
^ o P _o_ 

roc) 


- 


3h /en-K y[k, (Ah + gh /e )] 


o P 


Xh 


o 


O p 


r(K) 




+ k M 

V 


{l - » 


~Bh /e F[k, Xk h 


V o 


r(K) 


p 3k / e “1 

h t ® P 

y[K, (Xk^h^'+ 3h^/e^)]-j‘ 

L Xk h 

i* V o 

r(K) i 


Therefore, 


3E[E* ] 
A 
SM 


J(E, M, k^, h^) 


“• ’‘v’ ’'o' 


OM 


3E[E ] 



3M 


(5.12) 


(5.13) 


(5.14) 


(5.15) 
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For - 1, all M sensitivity comes from J(E, II, k^, > which will then 

have a mitiimum at M^, in which special case, according to Equations 
(5.11) and (5.12), 


3J(E, M, k , h ) 
’ v’ o 

3M 


at M = M 
o 


(5.16) 


As discussed in Section 4.4, evapotra aspiration efficiency, J, 

can be determined, for a given climate,' from observations of vegetation 

density and species by using the vegetal equilibrium hypothesis, Eq. 

(5.14). The actual procedure for doing this is to pick a value for the 

evaporation parameter, E, and calculate J for different values of M until 

evapo transpiration from soil moisture, Eq. (5.13), is minimized. If 

the vegetation density obtained which minimizes E„ is not equal to the 

observed value, E is incremented and a new M is found. For a fixed 

o 

climate, variations in E correspond to variations in the soil properties 
k-(l), c, n, and s^.' Therefore, what is actually done is seeking the soil 
which produces the observed vegetation canopy density for a specific 
climate. Once this value of evapotranspiration is found, the value of 
E is also knoxjn, which is a function only of the soil parameters for a 
given climate. 

With this information in mind, the following procedure is used 
to estimate the average areal effective soil parameters for a given set 
of climatic and vegetal parameters : 

1. A value for n is assumed and k is set = 1 

V 
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2 . 


3. 


4. 


5. 


. 6 . 


7a. 


7b. 


The above procedure is followed to determine E 

The lowest possible value for c, approximately 3.1, is 

picked as an initial value 

k(l) is calculated from Equation (5.3) 

With these values for the three soil parameters, n, k(l), 
and c, it can be seen from Eq. (4.21) that s^ remains as 
the only variable needed for determining E. With E known 
from step 2, is calculated 

Annual precipitation is calculated via Equations C5.5) 
through (5,10) 

If the annual precipitation from Step 6 is not equal to the 

actual mean rainfall, c is incremented upward from its 

initially low value and Steps 4-6 are repeated. 

Due to the approximation introduced by using Equation (5.3), 

the precipitation, P , calculated in Step 6 may never 

A 

exactly equal the actual mean value, m_ , for any value 

A 

as c is increased, coming to 
within AP^ of equality at intermediate c before diverging 
again for large c, Por low values of c, the calculated 
k(l) is large, representing a soil with high permeability 
and well connected pores. With evapotranspiration 
specified at the optimtim. (i.e, , minimum) yalue, a large 
precipitation is therefore calculated in order to produce 
the inevitably large groundwater yield of the highly porous. 


of c. P^ will approach 
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soil. For large c and small k(l), the soil is extremely 

impervious and the surface yield will be high. With 

minimum evapotranspiration, a large value for precipitation 

is again needed. Somewhere between these two extremes, a 

set of suitable soil parameters is obtained which gives an 

annual precipitation, which is closest to the actual 

mean, iil . This relationship is illustrated in Figure 5.4. 

Holding c constant at the value which gives the minimum 

AP., k(l) is then deviated from regression equation (5.3) 

A 

until another minimum in calculated precipitation is 
reached. If this value is above the mean precipitation, 
c is decreased, if it is below the mean, c is increased. 
Another search is done on k(l) until the minimum precipi- 
tation is found. This step is repeated until the minimum 
calculated precipitation is equal to the mean 
8. If the values obtained for k(l) and c are not consistent 
with the assumed porosity, n is adjusted to a more appro- 
priate value corresponding to a more pervious or impervious 
soil type depending on the values of k(l) and c. Steps 1 
through 7 are repeated. 

The soil parameters obtained from Steps 1-9 are used to 
construct the CDF of annual yield in the same manner as Eagleson (1978g) 
In this paper, the annual water balance is written as 
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Figure 5.4 

WATER BALANCE SOLUTIONS USING SOIL PROPERTIES FROM EQUATION (5.3) 





(5.17) 



E„ = R + R. 




In order to relate the annual water balance, (5.17), to the 

mean annual water balance, (5.1), Eagleson defines a climatic mean, m^, 

where = lOp and = E[E^ ], and expands (5.17) about this point in 
AAA 

a multidimensional Taylor expansion [Hildebrand, 1959, p. 353]. By 
taking expected values of this expansion term by term, neglecting higher , 
order terms, and assuming all variances, covariances, and curvatures 
are small, the "first order approximation" of E[Y^] is obtained: 


E[YJ - 



(5.18) 


This allows the use of the mean annual water balance equation to 
calculate annual values by letting the annual precipitation, and thus 
the average annual soil moisture, vary. The CDF of annual yield can 
then be calculated. Comparison of this CDF with that obtained from 
observations of annual streamflow provides the test for the accuracy 
of the estimated parameters, n, k(l) and c. 

Chapter 6 will present the results of this procedure in the 
form of annual CDF's of basin yield, in addition to verification of the 
equilibrium vegetation density hypothesis. 
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Chapter 6 


PRESENTATION OE RESULTS 

6.1 Verification of Vegetal Equilibrium Hypothesis 

The results of the applied methods of analysis explained in 
Chapter 5 are presented in this chapter- The vegetal equilibrium hypo- 
thesis is verified first in order to assure its validity for use in the 
estimation of soil parameters. 

Appendix A presents the individual catchments studied, the 
data used, location of the catchment, the values obtained for potential 
and actual evapo transpiration, vegetation density, and the estimated 
value of J- 

Eigure 6.1 presents the agreement of these experimental data 
with the hypothesized theoretical curves of Pigure 4.11. As can be 
seen, the dashed curves, which represent the derivation accounting for 
the presence of vegetation at the surface of the soil column in the 
exfiltration equation, provides a better fit for catchments with a 
vegetal canopy density greater than 0.2. This may mean that the presence 
of vegetation has a much greater effect on soil moisture exfiltration 
than previously believed. Although the equation used has serious flaws, 
they may be negligible compared to the possible importance of the 
presence of vegetation. 

Possible reasons for catchments W-4, W-5, and part of W-8 
lying above the curve may be ungaged yield which escapes through 
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Figure 6.1 

VERIFICATION OF VEGETAL EQUILIBRIUll HYPOTHESIS 
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groundwater aquifers, or flaws in the vegetal equilibrium hypothesis. 
Until it can be determined if all these yields are present in the 
observed streamflows, the vegetal equilibrium hypothesis would seem to 
give a reasonably accurate relationship between vegetation density and 
evapotranspiration . 

Also shoxm on Figure 6.1 are the results obtained from the 
empirical formulas developed by Czarnowski (1964) and Ritchie and 
Barnett (1971) . These functions exhibit the same type of relationship 
between evaporation efficiency and vegetation density, but do not fit 
the observed data quite so well as Equations (4.10) or (4.14). The fact 
that the data for these studies are primarily from agricultural areas, 
where cultivation and irrigation significantly violate the assumption 
of natural watersheds, is a likely reason for the poor observed fit. 

6.2 Estimation of Effective Average Areal Soil Properties 

To determine the accuracy of the procedure described in Section 
5.2, the two catchments studied here will be those studied by Eagleson 
(1978f, g); Clinton, Ma. and Santa Paula, Ca. Table 6.1 presents the 
list of necessary input variables (Eagleson, 1978g) and the computer 
program employed is listed in Appendix B. Tables 6.2 and 6.4 list the 
results obtained for the inputs given in Table 6-1. Listed probabilities 
are calculated for given P./nip using the Poisson model of Eagleson 
(1978b) . In Clinton, the value of is held constant for the entire 
range of soil moistures, while in Santa Paula, the vegetation density 
is allowed to vary with annual precipitation, as explained by Eagleson 
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(1978g) . These results are obtained using the form of surface runoff 
developed in Section 4.3 and presented in Section 5.2, (Equation 5.6). 
Figures 6.2 and 6.3 present the results in the fora of CDF’s of annual 
yield. 

Figure 6.2, ijhich represents Santa Paula, also shows the CDF 
obtained by Eagleson (1978g) for his silty-loam soil, which is listed in 
Table 6.3. The soil properties estimated by the algorithm explained in 
Section 5.2 indicated a slightly less permeable soil than the silty loam. 
This soil gives an improved fit over the entire range of CDF values, 
especially in the critical lower tail. 

The results for Clinton are illustrated in Figure 6.3. The 
soil properties obtained in this case indicate, again, a more impermeable 
soil than the silty-loam employed by Eagleson. Although these values 
for k(l) and c are quite different, the resulting CDF of annual yield is 
Indistinguishable from that obtained for silty-loam. To facilitate the 
comparison between the two results. Table 6.5 lists the annual water 
balance components for Clinton, using the silty-loam soil properties. 

Since the estimated soil properties represent a tighter soil which reduces 
the mobility of moisture, the soil moisture values are higher than for 
the silty-loam. The other major differences bet\-?een the two soils are 
the values for surface and groundwater runoff. The more permeable silty- 
loam yields a large groundwater component, and a surface runoff component 
which seems unrealistically low for all values of annual precipitation. 

In the case of the estimated soil properties, the surface runoff is 
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Clinton, Mass. 
.15 
3 

.32 

365 

.5 

.1 

0 

8.4 

94 

1 

.8 


cm/ day 

days 

'b 

days 

X 

days 

h. cm 
o 

w/e 

P 

T, °C 
A 

nip cm 
k 

V 

M 

o 


Santa Paula, Ca. 
.273 
10.4 
1.43 
212 


.25 

.1 


0 

13.8 

54 

1 


.4 


Table 6,1 

INPUT CLIMATE AND VEGETATION PARA^IETERS 
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EFFECriVli AREAL AVERAGE SOU PROPERTIES 


k(1)= 0, 2815708E-09 n=0* 35 c= 6.204 


so 

B 

Mo 

Pa/mpa 

Ya/mpa 

ESA 

EGA 

ETA 

PE OB 

0.30 

0.0 27 

0.068 

0. 2108 

0, 0158 

0.6264 

0.225 

10.533 

0.0216 

0.32 

0.038 

0.082 

0.2454 

0.0202 

0.7565 

0,335 

12. 158 

0.0314 

0.34 

0, 054 

0. 099 

Oc 2833 

0. 0258 

0.9071 

0.488 

13.901 

0,0446 

0.36 

0.074 

0.117 

0, 3246 

0,0329 

1.0807 

0. 696 

15.752 

0. 0619 

0.38 

0. 1 01 

0.138 

0.3694 

0.0417 

1.2798 

0.974 

17.694 

0,0840 

0.40 

0. 134 

0, 160 

0* 4177 

0, 0527 

1,5077 

1,339 

19.712 

0.1114 

0,42 

0.176 

0.183 

0,4697 

0, 066 3 

1.7679 

1.812 

21.786 

0.1446 

0.44 

0. 228 

0.209 

0. 5255 

0, 083 0 

2.0648 

2,413 

23.896 

0. 1839 

0,46 

0.293 

0.236 

0, 5854 

0. .1035 

2. 40 36 

3. 186 

26. 021 , 

0.2294 

0,48 

0,372 

0.264 

0,6496 

0.1285 

2.7909 

4, 148 

28, 141 

0,2809 

0, 50 

0. 468 

0. 294 

0-7187 

0.1589 

3.2345 

5.344 

30.233 

0,3380 

0,52 

0.582 

0.324 

0,7933 

0. .1956 

3.7443 

6, 816 

32. 277 

0.4001 

0.54 

0. 720 

0.355 

0,8741 

0,2397 

4.3323 

8.614 

34.253 

0. 4662 

0,56 

0, 882 

0. 367 

0,9621 

0. 2927 

5.0138 

10,794 

36.144 

0,5352 

0,53 

1 . 074 

0.419 

1 . 0585 

0. 3561 

5.8080 

13, 420 

37. 934 

0, 6056 

0,60 

1. 298 

0. 451 

1.1650 

0.431 5 

6.7388 

16.561 

39.610 

0.6753 

0,62 

1.560 

0.483 

1. 28 33 

0. 5210 

7, 8362 

20. 297 

41. 163 

0,7424 

0. 64 

1.864 

0.514 

1 ,4156 

■ 0.6269 

9.1381 

24.716 

4 2. 586 

0.8046 

0. 66 

2. 215 

0, 544 

1.5645 

0, 752 0 

10.6921 

29,914 

43.874 

0,8595 

0,68 

2,618 

0,574 

1.7331 

0, 8992 

12. 5580 

36.001 

45.027 

0.9054 

0.70 

3. 079 

0.602 

1 .9250 

1,0723 

14.81 16 

43.094 

46.046 

0,9412 

0.72 

3. 606 

0,6 29 

2t -1446 

1. 2754 

17.5485 

' 51.324 

46. 936 

0, 9668 

0.74 

4,204 

0,655 

2.3 968 

1,5134 

20.8905 

60.833 

47,703 

0.9833 


Table 6,2 

ANNUAL WATER BALANCE COMPONENTS, SANTA PAULA, CA, 
ESTIMATED SOIL PROPERTIES. EQUATION (5.6), = .4 



EFFECTIVE AREAL AVERAGE SOIL PROPERTIES 


Tt(1)= 0. 1200000E-08 n = 0.35 c= 6.000 


so 

E 

Ho 

Pa/mpa 

Ya/mpa 

RSA 

EGA 

ETA 

PROS 

0.26 

0.031 

0.073 

0. 2174 

0. 0106 

0.0538 

0.519 

11 .165 

' 0.0233 

0.28 

0.047 

0.091 

0. 2592 

0.0162 

Q, 0673 

0.809 

13. 120 

0.0360 

0.30 

0.068 

0.112 

0.3 062 

0.0242 

0.0835 

1 .224 

15. 225 

0-0539 

0.32 

0. 097 

0. 135 

0. 3586 

0, 0353 

0, 1027 

1,803 

17,456 

0.0784 

0.34 

0.135 

0.161 

0.4168 

0.0504 

0 .1256 

2.593 

19,789 

0.1109 

0.36 

0, 1 85 

0, 188 

0,4815 

0. 0705 

0,1527 

3 .654 

22,194 

0.1526 

0.38 

0,250 

0.218 

0. 5533 

0.0970 

0. 1850 

5,055 

24.639 

0,2046 

0. 40 

0.331 

0.250 

0.6332 

0.1315 

0.2234 

6.876 

27.09 2 

0. 2675 

0,42 

0. 433 

0. 284 

0. 7223 

0. 1756 

0.2691 

9.215 

29. 521 

0.341 0 

0.44 

0.559 

0.318 

0.8222 

0.2316 

0.3238 

12. 182 

31. 892 

0.4240 

0.46 

0. 714 

0.354 

0. 9347 

0.3018 

0.3896 

15.906 

34.177 

0,5142 

0,48 

0,902 

0.391 

1. 0620 

0, 3889 

0, 4688 

20. 533 

36,349 

0,6080 

0. 50 

1.129 

0.427 

1 .2070 

0.4962 

0,5647 

26.232 

38. 38 3 

0.7004 

0,52 

1.401 

0.464 

1.3729 

0. 6273 

0.6813 

33.191 

40.262 

0.7859 

0.54 

1.724 

0.500 

1.5634 

0.7061 

0. 8233 

41. 626 

41.973 

0. 8592 

0.56 

2. 1 06 

0.535 

1 .7830 

0.9773 

0.9968 

51 .777 

43.508 

0.9162 

0.58 

2. 554 

0.569 

2. 0 367 

1. 2059 

1.2092 

63.911 

44.864 

0.9559 

0.60 

3.078 

0.602 

2.3304 

' 1.4777 

1.4692 

78.327 

46, 043 

0.9799 


Table 6,3 

ANNUAL WATER BALANCE COMPONENTS, SANTA PAULA, CA. 
SILTY-LOAM SOIL PROPERTIES. EQUATION (5.6), M^ = .A 



ON 




EFFECTIVE AREAL AVERAGE SOIL PROPERTIES 


k(1)= 0, 1947604E-09 n=0.35 c= 7,399 


so 

E 

No 

Pa/mpa 

Ya/inpa 

ESA 

EGA 

ETA 

PROS 

0.60 

4. 387 

0. 800 

0. 6801 

0, ,1776 

7,5385 

9.153 

47.241 

0.0180 

0.62 

5.131 

0.800 

0.7225 

0.2145 

8.5004 

11, 666 

47.752 

0.0373 

0. 64 

6, 247 

0.800 

0.7725 

0.2598 

9,6665 

14,755 

48.195 

0.0766 

0.66 

7.560 

0,800 

0,8 317 

0.3151 

1 1, 0936 

18. 528 

48.562 

0, 1531 

0.68 

9. 097 

0.800 

0.9 023 

0.3826 

12.8580 

23,108 

48,852 

0. 2873 

0.70 

10.888 

0.800 

0. 9869 

0. 4649 

15,061 4 

28.636 

49.071 

0,4866 

0.72 

12.966 

0.800 

1,0 887 

0.5650 

17,8406 

35, 272 

49. 227 

0,7155 

0.74 

15, 367 

0, 800 

1.2118 

0,6870 

21 .3 8 03 

43 ,200 

49.332 

0. 8959 

0.76 

18. 130 

0.800 

1. 3612 

0, 8357 

25.9318 

52.623 

49.398 

0,9795 


Table 6.4 

ANNUAL WATER BALANCE COMPONENTS, CLINTON, MA. 
ESTIMATED SOIL PROPERTIES. EQUATION (5.6), = . 
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EFFECnVE aREAL AVERAGE SOIL PROPERTIES 


k(1)= 0. 1200000E-08 n=0, 35 c= 6.000 


SO 

E 

Mo 

Pa/mpa 

Ya/mpa 

RSA 

RGA 

ETA 

PHOB 

0,44 

5.540 

0.800 

0,7032 

0. 1933 

0, 2482 

17. 925 

47, 932 

0,0272 

0.46 

7. 074 

0.800 

0.7674 

0.2521 

0.2936 

23.404 

48. 441 

0,0716 

0.48 

8,940 

0.800 

0, 8446 

0. 3251 

0,3509 

30.213 

48.827 

0.1742 

0.50 

11.190 

0.800 

0,9 37 5 

0.4151 , 

0, 4239 

38, 598 

49. 099 

0.3673 

0,52 

13.884 

0.800 

1.0493 

0.5251 

0,5174 

48,839 

49. 27 4 

0. 6330 

0,54 

17.087 

0.800 

1, 1 837 

0. 6584 

0.6376 

61.250 

49,378 

0,8653 

0. 56 

20. 871 

0.800 

1 .3448 

0,8189 

0.7930 

76. 185 

49. 432 

0.9750 


Table 6.5 

ANNUAL WATER BALANCE COMPONENTS, CLINTON, MA, 
SILTY-LOAM SOIL PROPERTIES, EQUATION (5.6), = ,8 



% LESS THAN 



Figure 6.3 

FREQUENCY OF ANNUAL BASIN YIELD, CLINTON, MA 



greater and the groundwater runoff lower than for the silty-loam soil. 

The identical CDF’s of yield for the two soils can be explained in the 
following manner; 

In this development, the storage of moisture is not taken 
into account, therefore, yield is equal to precipitation minus evapotrans- 
piration. In the Clinton system, evapo transpiration is controlled primar- 
ily by the climate (Eagleson, 1978d) , and is relatively Insensitive to 
the soil properties except for extreme cases. Thus, for a given precipi- 
tation, evapotranspiration and hence yield, will be the same for different 
types of soil. The only variations occur in the proportioning of yield 
between surface and groundwater runoff. The permeable soil encourages 
gravitational percolation and hence groundwater, while the impermeable 
soil rejects precipitation as surface runoff. 

In Santa Paula, where evapotranspiration is primarily soil 
controlled, the yield is more sensitive to changes in the soil properties, 
and thus there is a difference in the CDF’s for the two different soils. 

In Figure 6.4, the estimated soil properties are used to show 
the effect on the yield CDF of the two methods of handling surface reten- 
tion in calculating surface runoff. As expected, the values obtained for- 
yield, using Eq. (5.5) are reduced from those calculated by Eq. (5.6) due 
to the reduction of rainfall excess in favor of surface retention. 

Although the difference between the two equations is not large. Equation 
(5.6) still fits the observed data better in the lower tail. 

Tables 6.6 and 6.7 list the CDF's obtained for Clinton and 
Santa Paula, using Eq. (5.5). Again, in the case of Clinton, the CDF is 
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EFPECTlVi: AREAL AVERAGE SOIL PROPERTIES 


k(1)= 0. 2964968E-09 n=0.35 6,423 


SO 

E 

«o 

Pa/mpa 

Ya/Vipa 

RSA 

EGA 

ETA 

PR OB 

0,32 

0.031 

0,074 

0, 2208 

0. 0128 

0,4166 

0,275 

1 1.'233 

0. 0243 

0,34 

0.044 

0,089 

0.2556 

0.0170 

0.5129 

0.40 6 

12. 885 

0.0348 

0.36 

0. 062 

0. 106 

0. 2937 

0.0225 

0.6267 

0.586 

14.649 

0,04S7 

0.38 

0,084 

0,125 

0,3353 

0.0295 

0. 7608 

0. 830 

16.513 

0, 0669 

0. 40 

0,112 

0.146 

0,3 803 

0,03 84 

0.9181 

1 ,153 

18. 463 

0,0898 

0.42 

0.149 

0, 168 

0. 4289 

0, 0496 

1,1018 

1.578 

20.481 

0.1182 

0.44 

0,1 94 

0,193 

0.4813 

0. 0638 

1.3161 

2. 127 

22. 548 

0,1525 

0. 46 

0,250 

0. 218 

0.5378 

0.0814 

1 .5658 

2.830 

24.646 

0. 1930 

0.48 

0,319 

0,2 46 

0.5937 

0. .10 33 

1,3567 

3. 720 

2 6,754 

0.2399 

0.50 

0.402 

0.274 

0,6645 

0.1302 

2.1960 

4.835 

28,8 49 

0,2930 

0. 52 

0.503 

0.304 

0.7357 

0, 1632 

2.5923 

6,221 

30.913 

0.3521 

0,54 

0,624 

0,334 

0.8131 

0. 20 34 

3.0565 

7.927 

32,924 

0,4165 

0.56 

0.768 

0.366 

0.8977 

0.2521 

3 ,6023 

10.013 

34.863 

0. 4852 

0,58 

0. 9 39 

0, 397 

0. 9908 

0, 3 1 0 9 

4,2466 

12,544 

36,715 

0, 5569 

0,60 

1,139 

0.429 

1 ,0939 

0,3816 

5.0109 

15.595 

38. 463 

0. 6297 

0. 62 

1. 37 4 

0. 461 

1.2087 

0.4662 

5,9222 

19,251 

40.097 

0,7014 

0,64 

1,647 

0,492 

1, 3375 

0, 5670 

7, 0 147 

23, 605 

41. 606 

0.7695 

0.66 

1.964 

0,523 

1 .4830 

' 0.6869 

8.3314 

28.763 

42.986 

0.8313 

0,68 

2. 329 

0.553 

1, 6482 

0. 8291 

9. 9270 

34.842 

44.231 

0.8842 

0,70 

2.748 

0,582 

1,8368 

0. 9971 

11, 8709 

41,972 

45, 342 

0.9266 

0.72 

3. 228 

0.610 

2.0531 

1.1953 . 

1 4.2509 

50.297 

46. 322 

0.9577 

0,74 

3,775 

0. 637 

2, 3024 

1. 4283 

17. 1789 

59,974 

4 7. 175 

0. 9783 


I Table 6,6 

ANNUAL WATER BALANCE COMPONENTS, SANTA PAULA, CA. 
ESTIMATED SOIL PROPERTIES, EQUATION (5.5), = ,4 



% LESS THAN 



Figure 6,4 

SENSITIVITY OF ANNUAL BASIN YIELD TO ITTO 1-lETHODS OF HANDLING SURFACE RETENTION, SANTA PAULA, CA, 



EFFECTIVE AREAL AVERAGE SOIL PROPERTIES 


k{1)= 0.1916398E-09 n=0,35 c= 7.664 


so 

E 

«o 

pa/mpa 

Ya/mpa 

ESA 

EGA 

ETA 

PROS 

0.64 

4.956 

0.800 

0.6876 

0.1805 

4.0637 

12.900 

47.669 

0,0207 

0.66 

6.023 

0,800 

0, 7389 

0. 2270 

5,0093 

16.331 

48.117 

0.0479 

0.68 

7.276 

0,800 

0, 8001 

0,2843 

6. 1903 

20. 530 

48. 493 

0. 1080 

0.70 

8,7 42 

0. 800 

0. 8736 

0,3545 

7,6819 

25 .637 

48.795 

0, 2275 

0,72 

10. 449 

0, 800 

0. 9620 

0.440 4 

9,5866 

31, 815 

49. 026 

0.4261 

0,74 

12. 428 

0,800 

1 .0690 

0.5457 

12.0451 

39. 249 

49. 194 

0. 6756 

0.76 

14.715 

0,800 

1 ,1990 

0.6745 

15. 2515 

48. 150 

49. 309 

0.8827 

0,78 

17. 345 

0, 800 

1.3576 

0,8322 

19.4752 

58,756 

49.383 

0,9786 


CO 

CjJ 


Table 6.7 

ANNUAL WATER BATJ^NCE COMPONENTS, CLINTON, MA. 
ESTIMATED SOIL PROPERTIES. EQUATION (5.5), = .8 



Identical to that obtained from Eq. (5,6). This can again be attributed 
to the fact that Clinton is primarily climate controlled, and evapotrans— 
piration is held almost constant near the potential regardless of the 
amount of water that is infiltrated or removed as surface runoff. 


In order to study the sensitivity of the results presented 
here to the vegetation density, values of that bracket the observed 
values are used in the soil property estimation program. Figure 6.5 
illustrates the results obtained for Santa Paula, which are listed in 
Table 6.8. Inputing an of 0.2 generates a set of soil properties 
that produces more yield and less evapotrans piration than the soil obtained 
using an of 0.4. By specifying such a low vegetation density, the 
vegetal equilibrium hypothesis used in the water balance produces a low 
value of evaporation efficiency, J (Figure 4.8). This corresponds to an 
annual evapotranspiration considerably below the potential. By reducing 
the evapotranspiration, the yield must be increased for a given precipi- 
tation, as can be seen by Eq. (5.1). 

On the other hand, attempting to input an which is greater 
than 0.41 does not give a solution. That is, no soil can be found for 
the Santa Paula climate which will produce a vegetation density much 


larger than the observed value of .4. 


The climatic variables. 



and m , at Santa Paula prohibit the system from sustaining a larger 
b 

vegetation density, and thus a higher evaporation efficiency. If annual 


precipitation is increased, or e^ decreased, the resulting increased 
availability of moisture would allow a greater M^. 

The same type of results are seen in Figure 6.6 and Table 6.9 


for Clinton. Even though the vegetation density is already large, and 
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EFFECTIVE AEEAI. AVERAGE SOU PROPERTIES 


lc(1)= 0. 9733S59E-09 n=0, 35 9.696 


so 

E 

Mo 

Pa/mpa 

Ya/mpa 

RSA 

RGA 

ETA 

PROB 

0.48 

0. 019 

0.056 

0. 2057 

0. 0349 

0.7780 

1 .105 

9.226 

0,0204 

0,50 

0, 026 

0.067 

0, 2401 

0.0476 

0. 9311 

1. 642 

10.390 

0. 0298 

0.52 

0. 034 

0,078 

0.2809 

0.0652 

1.1181 

2,401 

11.649 

0,0437 

0,54 

0. 045 

0.090 

0. 3298 

0. 089 1 

1. 3481 

3. 462 

13, 000 

0,0643 

0.56 

0. 059 

0.104 

0.3889 

0.1215 

1.6328 

4.926 

14.441 

0.0946 

0.58 

0. 077 

0. 120 

0.4607 

0, 1650 

1.9877 

6,923 

15.965 

0. 1386 

0,60 

0. 099 

0,136 

0.5485 

0. 2231 

2. 4329 

9.617 

17. 567 

0. 2010 

0.62 

0. 125 

0.154 

0.6564 

0.3002 

2.9938 

13.216 

19.238 

0. 2864 

0.64 

0. 158 

0. 174 

0. 7899 

0. 4016 

3.7033 

17.980 

20,968 

0.3972 

0,66 

0,199 

0.195 

0.9552 

0.5340 

4.6027 

24.231 

22, 747 

0. 5300 

0.68 

0. 247 

0. 217 

1. 1606 

0.7057 

5,7433 

32 .366 

24.563 

0,6726 

0,70 

0,3 06 

0.241 

1. 4159 

0. 9270 

7. .1878 

42. 870 

26,402 

0.8047 

0. 72 

0.376 

0.266 

1 .7333 

1 ,2101 

9.0110 

56.335 

28. 249 

0.9054 

0.74 

0, 460 

0, 292 

2, 1272 

1. 5699 

11.2994 

73,478 

30,091 

0.9652 


Table 6.8 

ANNUAL WATER BALANCE COMPONENTS, SANTA PAULA, CA. 
estimated soil properties, equation (5.6), = .2 
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Figure 6.5 

SENSITIVITY OF ANNUAL BASIN YIELD TO VEGETAL CANOPY DENSITY, SANTA PAULA, CA 







EFFECTIVE AFEAL AVERAGE SOIL PROPERTIES 


k{1)= 0.7135189B-09 n=0. 35 9. 247 


SO 

E 

Ho 

Pa/mpa 

Ya/inpa 

BSA 

EGA 

ETA 

PROS 

0,62 

1,700 

0,600 

0.7109 

0.2612 

6.8828 

17. 667 

42.272 

0.0309 

0,64 

2. 132 

0. 600 

0,8015 

0,3386 

8.1368 

23,695 

43,507 

0. 1097 

0.66 

2.655 

0,600 

0,9 140 

Oo 4387 

9. 7412 

31. 495 

44, 682 

0.3132 

0.68 

3.284 

0.600 

1.0539 

0.5672 

11.8042 

41,508 

45.758 

0,6433 

0,70 

4.037 

0.600 

1, 2281 

0, 7313 

14, 4697 

54.268 

46.703 

0,9109 


os 


Table 6.9 

ANNUAL WATER BALANCE COMPONENTS, CLINTON, 1-IA. 
ESTIMATED SOIL PROPERTIES. EQUATION (5.6), = .6 



EFFECTIVE AREAL AVERAGE SOIL PROPERTIES 


k{1)= 0, 7220880E-10 n=0.35 c= 5.449 


so 

E 

Mo 

Pa/mpa 

Va/mpa 

RSA 

EGA 

ETA 

PEOB 

0.56 

7,268 

0.900 

0,6513 

0.1321 

6.1125 

6,309 

48.805 

0.0103 

0.58 

8.7 31 

0.900 

0. 6741 

0. 1530 

6.7426 

7.639 

48.986 

0.0161 

0.60 

10,422 

0.900 

0.7000 

0. 1773 

7.4778 

9. 189 

49. 131 

0,0257 

0,62 

12. 370 

0. 900 

0.7295 

0,2056 

8.3421 

10.986 

49.242 

0,0415 

0.64 

14.602 

0.900 

0. 7633 

0, 2386 

9. 3664 

13, 061 

49.324 

0.0678 

0,66 

17. 1 48 

0,900 

0,8023 

0.2770 

10.5902 

15.446 

49. 382 

0. 1108 

0.68 

20. 043 

0.900 

0.8474 

0,3217 

12, 0647 

18.174 

49.420 

0. 1791 

0.70 

23,320 

0.900 

0.8998 

0. 37 38 ' 

13, 8567 

21,284 

49. 444 

0.2819 

0.72 

27, 018 

0.900 

0.9609 

0.4348 

16.0541 

24,816 

49. 458 

0,4235 

0,74 

31. 176 

0,900 

1. 0325 

0. 506 2 

. 18,7745 

28.812 

49.466 

0.5950 

0.76 

35.837 

0.900 

1.1167 

0.5904 

22. 1777 

33. 318 

49, 470 

0. 7672 

0.78 

41.046 

0.900 

1.2164 

0.6901 

26.4843 

38.384 

49.472 

0,9003 

0,80 

46. 851 

0.900 

1.3355 

0. 8092 

32. 0064 

44. 062 

49. 472 

0,9721 


Table 6.10 

ANNUAL WATER BALANCE COMPONENTS, CLINTON, MA. 
ESTIMATED SOIL PROPERTIES. EQUATION (5.6), = . 
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evapo transpiration is near the potential, it is still impossible to find 

a soil x^hich alloxirs an M much larger than the observed value of 0.8. 

o 

Again, reduction of produces a soil which generates a larger amount 
of annual yield for the same reasons mentioned for Santa Paula.. 

On the basis of these comparisons xje see the soil properties 
determined from the estimation algorithm describe the behavior of these 
two systems very well through the water balance model. A brief summary, 
and conclusions draxm from these results will be presented in Chapter 7. 

Although the yield CDF's for Clinton derived from varying 
soil properties are identical, the values obtained for the average annual 
soil moisture vary significantly between the silty-loam soil and the 
soil found from the algorithm. Since soil moisture is a state variable, 
it is desirable to be able to verify the accuracy of its prediction. 

One possible method for doing this would he to compare the CDF's of 
surface runoff, rather than total yield. It has been noticed that the 
surface runoff components of the annual water balance are much more sensi- 
tive to changes in soil properties than is the total yield. One problem 
with this, hox^ever, is the lack of measurements of surface runoff, al- 
though streamflow in arid climates may actually be composed totally of 
surface runoff. 
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Figure 6.6 

SENSITIVITY OP ANNUAL BASIN YIELD TO VEGETAL CANOPY DENSITY, CLINTON, MA 







Chapter 7 


SUMMARY AND CONCLUSIONS 


A one- dimensional water balance model (Eagleson, 1978a, b,c,d, 
e,f) is employed to parameterize the climate-soil-vegetation relationship 
at the land-air interface. A vegetal equilibrium hypothesis proposed by 
Eagleson (1978f) provides a second relationship between the climate, soil 
and vegetation. 

Improvements are made in the method of calculating the bare 
soil component of evaporation, and in the way surface retention is 
handled. 


The vegetal equilibrium hypothesis is developed, and its use 
in the water balance is explained. The sensitivity of this hypothesis 
to various parameters of the evapotranspiratlon function is explored. 

It is found that the two parameters to which the system is most sensitive 
are Bh /e , which can be readily evaluated, and k , whose value is uncer- 
tain. It is believed that is usually equal to one, except in very 
dry climates, where the plants transpire at a rate less than an equiva- 
lent area of bare wet soil. In this work, k is held at its nominal 

V 

value, 1. 


Reasonable verification of the vegetal equilibrium hypothesis 
is obtained through comparisons of the theoretical relationship between 
density of canopy cover and the evapotranspiratlon efficiency to data 
obtained from observations in watersheds representing various types of 
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clima-tes . 


An algorithm is derived which searches for the soil properties 
that produce, in a given climate, the level of evapotranspiration deter- 
mined through observations of vegetation density. By using the vegetal 
equilibrium hypothesis, the water balance, and a regression equation relat- 
ing the soil’s intrinsic permeability and pore size distribution index, 
a consistent set of soil properties is found which generates the implied 
evapotranspiration and also satisfies the mean annual water balance. 

This estimation of soil properties produces results, through 
the water balance, in the form of CDF’s' of annual basin yield, that 
describe the observed behavior of the Clinton and Santa Paula systems 
very well. In both Clinton and Santa Paula, the soils determined were 
slightly less permeable than the silty-loam which Eagleson (1978g) used 
as his best-fitting soil. These soils also produce a more realistic 
(although unverified) surface runoff component than those used by 
Eagleson. 

A remaining important question is the sensitivity of the 
water balance model to the vegetation parameters, and k^. Inclusion 
of this analysis was beyond the scope of this study, and it is left as 
an important subject of future work. 

From this summary, the following conclusions may be dra^TOi: 

1. The vegetal equilibrium hypothesis is sufficiently valid 
to justify its use as a supplementary water balance rela- 
tionship between the soil, climate, and the vegetation. 
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2. The algorithm for estimating the effective areal soil 
properties works well, producing CDF's of annual yield 
which fit the observed CDF's closely. 

3. It is more accurate to subtract surface retention from 
the volume of infiltrated precipitation at the beginning 
of the rainfall period than from the rainfall excess. 

4. Use of the vegetal equilibrium hypothesis and the soil 
estimation algorithm should facilitate the incorporation 
of the areal variability of soil properties and soil 
moisture dynamics into global climate models . 
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Chapter 8 


RECOMMENDATIONS FOR FUTURE WORK 


Questions remaining and subjects for future study are: 

1. Evaluation of the Philip exfiltration equation for a 
varying soil column cross-section. 

2. Sensitivity of the water balance to vegetation through 
the parameters, and k^. 

3. Development of a procedure for determining the accuracy 
of predicted values of average annual soil moisture. 
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Appendix A 


DATA FOR CATCHMENTS STUDIED IN VERIFICATION OF 
VEGETAL EQUILIBRIUM HYPOTHESIS 


98 



W-1 


Location: Albuquerque, New Mexico 

Latitude: 4* “ 35 °N 

Rainfall; P = 4.37 in. 

Streamflow: Q = .2 in. 

Season Length = 4 mos., July — Oct. 

Cloud Cover: N = .37 

Humidity; S = 39.97% 

Temperature; T = 69.6l°F 

Vegetation Density; M = .12 to .15 

o 

Albedo: A = ,25 to ,3 

e = 15.47 in/season to 14.19 in/season 

P 

J = .27 to .294 

Watershed Conditions; Rough broken rangeland. About 85% is bare. Sparse 
vegetation consists of short grasses, shrubs, and a few small juniper 
and pinion trees. 

Comments: The value for is estimated directly from the percent bare 

ground, and taking into account the cro\m spread of the trees. 

Source*: Hydrologic Data for Experimental Watersheds in the United States, 

1967. U.S.D.A. 

* Indicates reference from which vegetation density values are obtained, 
and in some cases, precipitation and streamflow data as well. All other 
data is obtained from U.S. Weather Bureau publications and U.S.G.S. 
reports of surface water resources. 


E = 4.37 - .2 
A 

= 4.17 in. 
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W-2 

Location: Cornfield Wash, New Mexico 

Latitude; ^ = 35°N 

Rainfall: P =» 6.29 in. 

Streamflow: Q = function of M 

o 

Season Length = 4 mos-, July - Oct. 
Cloud Cover; N = .37 

Humidity; S = 39.97% 

Temperature; T = 69.61°F 


Vegetation Density: 

M 

o 

= .16 

Q 

=1.07 in. 


= 5.22 in. 


M 

o 

= .24 

Q 

= .28 in. 


= 6.01 in. 

Albedo: 

A 


.25 

to 

.30 



e 

P 

= 

15.47 

in/ season to 

14.19 

in/ season 

M = .16 
o 

J 

= 

.34 

’ to 

.37 


M = .24 
o 

J 


.39 

to 

.42 



Watershed Conditions; The dominant vegetation is galleta grass. Remain- 
ing areas have a mixture of other grasses, Russian thistle, and big sage- 
brush in small upland drainages. 

Comments: Runoff data was recorded as a function of percent hare soil in 

the paper used as the source, therefore, the calculation of E gives two 
values, one for each and Q data pair. Vegetation density values were 
recorded for each value of percent hare soil, and the two extreme values 
were used here. 
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W-2 (continued^ 


Source: 


F. A. Branson, and J. B. Owen» "Plant Cover, Runoff and Sediment 
Yield Relationships on Mancos Shale in Western Colorado," W.R.R. , 
6 ( 3 ), 1979 . 
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W~3, W-4, W-5 


Location: 

latitude: 

Rainfall: 

Streamflow: 


Season Length 
Cloud Cover: 
Humidity: 
Temperature: 


Tombstone, Ariz. 

(j> = 32'‘N 

= 8.01 in. 
= 8.44 in. 
= 7.56 in. 

= 3 mos., July ~ Sept. 

N = .35 

S - .4467 

T = 82.17°F 


P ^_3 == 8.65 in 
^W-3 ” *64 in. 

V4 " 

%-4 ~ -21 in. 

P-, , = 8.65 in. 
W-5 

= 1.09 in. 




V: 




Vf 


Vegetation Density: M 


M 

.35 

to 

.4 

°W-3 




M 

.25 

to 

.3 

°W-4 




M 

.2 

to 

.25 

W-5 





Albedo: 


A = 

.24 

to 

.30 

e = 
P 

13.45 in/season 

to 

12.14 

“^W-3 "" 

.60 

to 

.66 

■^W-4 

.63 

to 

.70 


.56 

to 

.62 


Watershed Conditions: All watersheds have cover of desert shrubs 

(whitehorn, creosote bush, tarbush) with an understory of grass (black 
grama, tobosa grass, blue grama, sideoats grama, and curly mesquite grass). 
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W-3: Entire area covered by shrubs with 38% crown spread. = .35 to .40. 

W-4: 78% of area covered by shrubs with crown spread of 30%. Remaining 

22% covered with grass with .2% basal area. ~ ,25 to .3. 

T^-5: Shrub canopy approximately 20%, Remaining area covered by grass 

with .2% basal area. M ~ .2 to .25. 

o 

Comments: The three watersheds are all sub-catchments of a larger catch- 

ment. Therefore, while vegetation densities and streamflow vary slightly, 
the annual climatic properties are all the same. 
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W-6 

Location: Flagstaff, Arizona 

Latitude <|) = 35 “N 

Rainfall: P = 12.38 in. E = 11.98 in. 

A 

Streamflow: Q » .4 in. 

Season Length = 7 mos., July - Jan. 

Cloud Cover: N .4 

Humidity: S = .52 

Temperature: X == 47.13°F 

Vegetation Density: = .3 to .35 

Albedo: A = .2 to .25 

e = 20,63 in/season to 18,95 in/season 

P 

J = .58 to .63 

Watershed Conditions: The terrain is undulating uplands dissected by many 

small drainages. The vegetation is mainly upper pinion juniper woodland 
with a sparse understory of grasses. 

Source: Bro^m, H.W., "Characteristics of Recession Flows from Small Water- 

sheds in a Seraiarid Region," W.R.R. , 1(4), 1965. 
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W-7 


Location: 
Latitude: 
Rainfall: 
Streamflow: 
Season Length 
Cloud Cover: 
Humidity; 
Temperature; 


Badger Wash, Colorado 
(j) = 38“N 
P - 4.69 in. 

Q = function of M 

= 6 mos - , August . - Jan . 
N = .5 

S = .4817 

T = 47.8“F 


Vegetation Density; 

Mo = .13, 

Q = .96 in. , 


E = 3.73 in 
A 


M ^ .26, 

Q = .35 in.. 


E - 4.34 in 
•^A 

Albedo: 

A = 

.25 

to 

.30 


e = 
P 

16.04 in/season 

to 

14.66 in/season 

- 

M = .13: 
o 

J = .23 . 

to 

.25 


M = .26: 

J = .27 

to 

.30 


Watershed Conditions: The catchment is in a semiarid area with pre- 

dominantly desert-type shrubs. 


Comments: This data was obtained in the same way as that for W-2. Thus, 

the values for J are presented in the same way. 

This watershed is located in an area where there is considerable 
snowfall. The model used in this work does not account for snowmelt in any 
way, and only works with yield resulting from precipitation in the form of 
rainfall. Therefore, if the yield measurement includes runoff from snow- 
melt, the value of precipitation used here is not large enough to account 
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W-7 (continued) 

for that much streamflow, and the resulting calculated value of actual 
evapotranspiration is too small. It would not be surprising then if the 
value plotted for J vs, is below the hypothesized curve. 

Source; Branson, F, A. and J. B, Owen, "Plant Cover, Runoff, and Sediment 
Yield Relationships on Mancos Shale in Western Colorado," W.R.R. , 6(3), 
1979. 
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W-8 


Location: 

Santa Paula, 

California 


Latitude : 

= 34.4°N 



Rainfall: 

Streamflow: 

P = 21.26 in 
Q = 6.85 in. 

• 

= 14.41 in. 

Season Length: 

= 7 mos . , 

Oct. - Apr. 


Cloud Cover: 

N = .37 



Humidity: 

S = .6897 



Temperature: 

T = 53.06'’F 



Vegetation Density: = 

.35 

to .5 

Albedo: 

A - 

.2 

to .32 


e = 

P 

21.23 in/season 

to 16,73 in/season 


J = 

.68 

to .86 


Watershed Conditions: Fairly rugged terrain with wide variation of vege- 

tation type. Dominant species are desert-type shrubs which are common in 
Southern California mountain ranges. 

Source: 1) Eagleson, P. S., "Climate, Soil and Vegetation," Parts 1-7, 

W.R.R. , 14(5), Oct. 1978. 

2) On-site observations. 
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W~9, W-10 


Location: 

Chickasha, Oklahoma 



Latitude : 

(j) = 




Rainfall: 

Streamflow; 

P n - 23.52 in. 
w-9 

- 1.12 in. 


%~9 

22.40 in. 


P-T in 23.52 in. 
W“1U 

Vio - 3-” 


\ 

^-10 

= 19.75 in. 

Season Length: 

~ 7 mos., Apr. - 

Oct. 



Cloud Cover: 

N = .47 




Humidity: 

S = 67% 




Temperature: 

T = 70.61“P 




Vegetation Density: M = 

°W-9 

M 

°W-10 

.45 

-2 

to 

to 

.57 

.3 

Albedo: 

A 

.18 

to 

.24 


e 5= 

P 

28-80 

in /season to 

26.09 in/season 

W-9: 

J 

.78 

to 

.86 

W-10: 

J 

.69 

to 

-76 


Watershed Conditions: The vegetation of both catchments consists of native 

grasses (buffalo grass, blue grama, little bluestem). Values of II are 

o 

interpreted from radiation shielding values obtained from average values of 


leaf area index and percent mulch cover. The equation used is (2) 


R -.4(L,,+2.5M) 
~ , ns Ax 

M 1 — := — = e 
o R 

no 
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where 


R = net radiation reaching the soil surface 
ns 

R = net radiation above plant canopy 
no 

L. . = leaf area index 
Ai 

M = fraction of surface covered by Mulch 

Comments: It is reported in the source paper that ¥—10 is constantly 

overgrazed, thus, it is likely that the value obtained for is unnaturally 
small, and the plotted position of this catchment will be above the hypo- 
thesized curve. 

Source: 1) Hydrologic Data for Experimental Agricultural Watersheds in the 

U.S. 1976. U.S.D.A. 

2) J. T. Ritchie, E. D. Rhoades and C. W. Richardson, "Calculating 
Evaporation from Native Grassland Watersheds," Transactions of 
the A.S.C.E ., Aug. 1976. 
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W-11 


Location: Clinton, Massachusetts 

Latitude: ^ = 42.50°N 


Rainfall: 

P = 43.82 in. 

= 22.01 in. 



^A 

Streamflow: 

Q - 21.81 in. 


Season Length: 

= 12 mos. 


Cloud Cover: 

N = .35 


Humidity: 

S = .70 


Temperature: 

T = 47.12“? 


Vegetation Density: ~ 

to . 9 

Albedo: 

A = .25 

to .30 


e = 24.25 in/season 

P 

to 21.64 in/season 


J = .91 

to 1.02 


Watershed Conditions: No specific conditions are available, only the range 

of vegetation density. 

Source: 1) Eagleson,,P. S\ , ""Climate, Soil and Vegetation," Parts 1-7, < ■ 

' 14(5), Oct. 1975. 

2) Visual observations of nearby "watersheds. 
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Appendix B 


FORTRAN PROGRAM FOR ESTIMATION OF SOIL PROPERTIES 
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c THIS PROGRAM CALCULATES EFFECTIVE AREAL AVERAGE SOU PROPERTIES. VJHEN 
c THE SOIL PROPERTIES ARE VARIED USING A REGRESSION EQUATION, CALCULATED 
e PRECIPITATION, Pa, REACHES A MINIMUM AT A MEDIAN VALUE OF c. THE PARA- 
c I'ETER, k(1), IS THEN DEVIATED FRCM THE REGRESSION UNTIL ANOTHER MINI- 
c MUM Pa IS FOUND. DEPENDING ON WHETHER THIS VALUE FOR Pa IS ABOVE OR 
c BELOW THE KNOWN VALUE OF mpa, THE PARAMETER c, IS INCREMENTED DP OR 
0 DOra, AND k(1) IS SEARCHED AGAIN UNTIL ANOTHER MINIMUM IS REACHED. , 
e THIS INCREMENTATION AND SEARCHING IS CONTINUED UNTIL THE MIMIMUN Pa * 
c FOUND IS EQUAL TO mpa. 


integer change , f tm , cf bl , runs, number , mon, it er 

real*8 ranu,p1 

real mtb,mtr,mh,inpa 

real mi,iao,m,n,nu,k1 ,k2 

0 DIMENSIONLESS INFILTRATION DIFFUSIVITI 

fii(d,so)=1,/(d*(1.-so)»»(1.45-.0375*d)+5./3.) 

C PORE SHAPE PARAMETER 

fi(em)=10.»*(.66+.55/em+.l4/em*»2.) 


print, 'Input parameters in the correct units. * 
print, *ep,em/day uitb,days mtr,days tau,days kappa,-.* 
print, *ho, cm w/ep,- ta, degrees C. * 
input , epr , mtb, mtr , t au , ak , ho , wep , ta 

pistol=1 

c IF pistols 1, THE ARRAY OF FACTORIilLS IN THE CDF SUBROUTINE HAS NOT 
c BEEN CALCULATED YET. ONCE pistols2, THE FACTORIALS HAVE BEEN STORED 
C AND THE LINES WHICH DO THIS CALCULATION ARE THEN SKIPPED. 

5 print, 'Input mpa, cm kv,- Ho,- n,- 

input , mpa , akv , mo , n 

print, 'For annually varying Mo, type 1, for constant Mo, type 2* 

input, mon 

if(mon.eq.O)stop 
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changes 1 

c IF changes 1, SOIL PROPERTIES ARE ROT YET DETERMINED. IF change=2, THE 
e SOIL PROPERTIES HAVE BEEN DETERMINED AND ONLY THOSE STEPS NEEDED FOR 
c DERIVING THE CDF OF THE WATER BALANCE COMPONENTS ARE USED. 


runs=1 j 
G IF runs=1, THIS IS THE FIRST SET OF SOIL PROPERTIES USED, AND NO COM- 
c PARISON OF CALCULATED Pa IS POSSIBLE. IF runs=2, THE NH^ VALUE OF Pa 
c IS COMPARED TO THE OLD VALUE TO SEE IF A MINIMUM HAS BEEN REACHED. 


cfbl=1 

c IF cfblsl, THIS IS THE FIRST DEVIATION OF k(1) FRCM THE REGRESSION 
G AND NO COMPARISON OF Pa IS DONE. IF cfbl=2, THE Pa CALULATED WITH 
e THIS k(1) IS COMPARED TO THAT CALCULATED USING THE PREVIOUS k(1) TO 
e SEE IF THE SECOND MINIMUN HAS BEEN REACHED. 


ftm=1 

e IF THE SECOND MINIMUM HAS JUST BEEN FOUND, BUT IF THIS MINIMUM 

.c Pa?Jmpa, e MUST BE CHANGED AND THE ENTIRE PRWESS MUST BE REPEATED, 
c THE VALUE OF THE DIFFERENCE BETWEEN THE MINIMUM Pa AND mpa, awbal, IS 
e PRESERVED AND COMPARED TO THE NEXT ONE OBTAINED. THIS COMPARISON IS 
e SIGNALED VJHEN ftra=2. WHEN awbal < .001, THE SOIL PROPERTIES HAVE BEEN 
G FOUND. 

c SET INITIAL VALUES 
Pl=0.0 
so=0.0 
dcs= . 1 
dics=.1 


c COMPUTE WATER CONSTANTS 
c sut=SURFACE TENSION 

c nu=VISCOSITI 

c gamsw=SPECIFIC WEIGHT 
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call VJATCN(ta,sut,nu,gainsw) 

e COMPUTE CLIMATIC PARAMETERS 

delta=1 ./mtr 
mh=nipa/ ( tau/ ( mtb+nitr ) ) 
mnu=tau/( rafcb+mtr ) 
mi=mh/mtr 
eta=1 ,/mh 
alpha=1 ,/mi 
pi=3. 14159 
beta=1 ./mtb 

epa=epr*tau*mtb/ ( mtb+mtr) 
al=ak/ioh 
alhsal*ho 
bhe=beta*ho/epr 
if Cho,eq.0.0)goto 10 
ble=bet^(al*epr) 
goto 20 
10 ble=0.0 

20 alkh=alh*akv 

blke=ble/akv 

30 if ( change. eq . 1 ) goto 40 

print 

print 

print,* so E Mo Pa/mpa Ya/mpa RSA 

do 400 i=1,45 

so=so+,02 

e= e cn st *so ** d2 

fiid=fii(di,so) 
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40 


if (change. eq.2)goto 45 
goto 50 
45 . if (mon.eq.2)goto 60 

goto 55 

c TO SPEED UP THE SEARCH FOR THE VALUE OF e THAT MINIMIZES ETA AT THE 
c OBSERVED Mo, e AND m ARE GIVEN INITIAL VALUES DEPENDING ON THE VAL- 
G UE OF Mo. BY PICKING A VALUE FOR e, THE m THAT MINIMIZES ETA CAN BE 
o FOUND. IF THIS m;?Mo, ANOTHER e IS PICKED UNTIL m=Mo. 

50 if(mo.ge..2)e=.3 

if(mo.ge,.3)e=,5 
if(mo,ge, ,4)e=1 . 
if(ino.ge..6)e= 3. 
if(mo.ge..7)e=6. 
if(mo.ge..8)e=10. 
if (mo.ge. . 9)e=20 , 

55 if(e.ge..01)bm=.1 

if(e.ge..1)bm=.4 
if(e.ge.1.)bra=.6 

if (e.ge. 10,)bm=,9 
if (mo.lt,.4)de=.01 
if (mo.ge.. 4)de=.1 
number=1 
60 iter=1 

dm= , 0 1 

70 bmkv=bm*akv 

c COMPOTE EVAPOTRANSPIRATION PARAMETERS, B & C. 

. b=((1 .-bm)/(1 .+bmkv)+(bmkv*bm)/(2,*(1 .+bmkv)**2.)) 

if (bmkv.eq.O.O)goto 80 
c=1 ./(2.*(brakv*bmkv)) 
goto 90 
80 c=1.e10 

90 be=b«e 
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ce=aminl ( 0 * 6 , 80 .) 

gamk=gafflt(ak, alh)/gamnia(ak) 
gamkl=gamt ( ak, alh+bhe )/ gamma ( ak) 
gambe=gamt( 1 .5 ,be) 
gamce=gamt ( 1 . 5 , ce ) 
gamkv=gamt ( ak, alkh) /garama(ak) 
gamkvl=gamt ( ak, ( alkh+bhe) ) /gamma ( ak) 

cSssftswsssf-a^sssssssasstKSS**®****^***********^****^?*************#******** 
c COMPUTE ANNUAL EVAPORATION FROM SURFACE RETENTION 

era=epr/beta*((1 .-bm)*(1 .•-exp(>-bhe)*(1 .-gamk)-(1 ,+ble)**(-ak) 

& *gamkl) 

& +bmkv*(1 ,-exp(-bhe)*(1 .-gamkv)~n .+ble)**(-ak)*gamkvl))*mnu 

eram=era 

c COMPUTE INTERSTORM BARE SOIL EVAPORATION 

es j=gamk“( 1 ,+ble) **( -ak) *gamkl*exp ( -be)+ 

& (1 .-gamk) * ( 1 . -exp ( -be-bhe ) * ( 1 .+bmkv+sqrt ( 2 . *b ) *e-wep ) 

& +exp(-ce-bhe)*(bmkv+sqrt(2.*c)*e) 

& +sqrt(2.*e) *exp ( -bhe ) * ( gam ee-gambe)) 

& +(1 .+ble)**(-ak)*gamkl*(sqrt(2.*e)*(gamce-gambe) 

& +exp(-ce)*(bmkv+sqrt(2.*c)*e) 

& -exp(-be)*( bmkv+s qrt ( 2 . *b ) * e-wep ) ) 


0 COMPUTE EVAPOTRANSPIRATION FUNCTION 

hj=1 ./(1-bm+bmkv)*( (1-bm) *esj+bmkv) 
ETN =h j * ( 1 . -bm+bmkv) 


if (change. eq. 2) goto 95 
goto 100 

95 if(mon.eq.2)goto 160 
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e THESE LINES FIND THE M THAT MINIMIZES ETA. 

c IF iter=1 , IT IS THE FIRST TIME THROUGH AND NO COMPARISON IS MADE 
100 if (iter.eq.l)goto 120 

if (absCdm). It.. 000001) goto 150 
if(ETN.gt.ETMIN)goto 110 
goto 120 

110 bm=bm- 1 . 5 dm 

dm=-.5*dm 
goto 130 
1 20 ETMIN=ETN 

bmin=bm 
iter=2 
bm=bm+dm 

130 if(bm)1il0,70,70 

1^10 bm=.1*(bm-dm) 

q=q+1 

if(q.lt.4)goto 70 

c AT THIS POINT, NO Mo CAN BE FOUND THAT IS GREATER THAN 0,AND NEl^ PAR- 
c AMETERS MUST BE INPUT. 

goto 395 

150 bmsbmin. 

ETN=ETMIN 

160 if (change. eq.2)goto 230 

c THESE LINES FIND THE E CORRESPONDING TO THE GIVEN Mo. 
c IF nuraber=1, IT IS THE FIRST TIME THROUGH AND NO COMPARISON IS MADE, 
diffsmo-bm 

if(abs(diff).lt..0001)goto 200 
if(number.eq.1)goto 170 
if(diff»diffold.le.O.O)goto 190 
if (number, eq. 2) goto 180 
170 if(diff.lt.0.0)de=-1.»de 
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number=2 

180 diffold=diff 
e=e+de 
goto 60 

190 de=-de».5 

diffold=diff 
e=e+de 
goto 60 

200 continue 

c AT THIS POINT, THE VALUE OF e HAS BEEN DETERMINED, AND SOIL PROPER- 
c TIES ARE NOW SEARCHED. 

cs=4. 

210 m=2./(cs-3.) 

fie=fi(in) 
dE=2.+1 ./m 
di=cs-1 ./m-1 . 
d2=dE+2. 
fied=fie(dE) 

c REGRESSION EQUATION 

k1=(m/512.7)*»2,75 

k2=k1 

dk1=k1/10. 

220 continue 

bkl =k1*gamsw/nu 

sil =sqrt (n/( k1 *f ic) ) *sut/gamsw 
sigc=n*eta®*2 . *bk1 *si1/ ( pi-m*delta) *720 00 . 
ecnst=2 . *beta*n*bk1 *si1 ied/ ( pi*m*epr**2 . ) *86400 . 

c SOIL MOISTURE IS CALCULATED. 

so=(e/ecnst)**(1 ./d2) 
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fiid=fii(di,so) 


o COMPUTE WATER BALANCE 

o COMPUTE ANNUAL E?AP0 TRANSPIRATION 

230 ETA=ETN*epa 

c**»»»»*a»»s»»»»sa»9s****»»»«**»«»*s**»*a»*s**»«*a****»»*»»»*»****»»»«»s 
c COMPUTE ANNUAL GROUNDWATER RUNOFF 

RG A= tau*bk 1 *so cs *86 40 0 


sigrf = (sigc*fiid*( 1 ,-so) **2. ) **.33333 
g=alpha*bk1 *86400*. 5*( 1 .+so**cs) 
blop=g+2 , *sigrf 
if(blop.gt.85.)blop=85. 

blip=exp(-blop)*gamma(sigrf+1 .)*sigrf**( -sigrf) 
if(blip.gt..95)blip=.95 

c COMPUTE PRECIPITATION, YIELD, RUNOFF 

Pa=(ETA+RGA)/(l .-blip) 

RSA=blip*Pa 

Ya=RSA+RGA 

if (change. eq.2)goto 380 
awbal=Pa-mpa 

G NOTE-awbal IS THE DIFFERENCE BETOEEN CALCULATED Pa AND KNOira mpa. THE 
e FOLLOWING LINES WILL PERFORM THE SEARCH FOR SOIL PROPERTIES THAT PRO- 
c DUCE Pa-mpa. 

if (cfbl.eq.2)goto 260 
if (ftm,eq.2)goto 280 
if (runs. eq.1) goto 250 

C THESE LINES PERFORM THE FIRST MINIMIZATION WHICH ADHERES TO THE RE- 
c GRESSION EQUATION, 
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if Cabs(dcs) .It. .OODgoto 260 
If (awbal.gt.awbol)goto 240 
goto 250 

240 cs=es-1 .5®dcs 

dcs=-.5*dcs 
goto 210 

250 awbol=awbal 

cs=cs+dcs 
runs=:2 
goto 210 

C THESE LINES PERFORM THE SECOND MINIMIZATION WHICH HOLDS C CONSTANT 
c AND DEVIATES k(1) FROM THE REGRESSION. 

260 if(cfbl.eq.2)goto 2?0 

if(abs(awbal).lt..001)goto 320 
if(cfbl.eq.1)goto 280 

270 if(abs(dk1).lt.k2/1000.)goto 320 

if (awbal.gt.awbol)goto 290 

280 awbol=awbal 

k1=k1-dk1 

c SINCE k(1) VARIES BY ORDERS OF MAGNITUDE, dkl MUST BE REDEFINED IF 
c k(1) GETS TOO BIG OR SMALL. 

if(k1.1t.k2/9)goto 300 
if (k1 .gt.k2*9)goto 310 
cfbl=2 
goto 220 

290 k1=k1+1,5*dk1 

dk1=-.5*dk1 
goto 220 

300 dk1=dk1/10. 

k2=k1 
goto 220 
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310 dk1=dk1»10 

k2=k1 
goto 220 

e THESE LINES PERFORM THE SEARCH ON c MICH LOCATES THE MINIMUN Pa FROM 

0 THE ABOVE PROCEDURE WHICH EQUALS mpa. 

320 efbl=1 

if(abs(awfaal),lt..01)goto 360 
if (ftm.eq.l)goto 330 
if (awbal®awbold.lt.0.0)goto 350 
goto 3^0 

330 if (awbai.gt.0.0)dics=-1 .*dic5 

ftm=2 

340 awbold=awbal 
cs=cs+dies 
goto 210 

350 dics=“dics®.5 
awbold=awbal 
G3=cs+dics 
goto 210 

360 print, ’AVERAGE EFFECTIVE PARAMETERS’ 
print 370,e,so,k1 ,es 

370 format( ’E=' ,f6.3,2x, *so=’ ,f5.3,2x, ’k(1)=’ ,el6.7,2x, ’e=» ,f6.3) 
ehange=2 
so=0.0 
goto 30 

380 y1 =Ya/mpa 
p1=Pa/mpa 

if(p1.1t..2)goto 395 

e COMPUTE CDF OF PRECIPITATION 

call PROBZ(mnu,p1 ,prob,ak) 
if(prob.lt..009)goto 395 


120 



if Cprob.gt. ,99)goto 410 
print 390 ,so,e,bm,p1 ,y1 ,RSA,RGA,ETA,prob 
390 format(f4.2,3x,f6.3,3x,f5.3,3x,fT.4,3x,f7.4,3x,f7-4,3x,f7.3,3x 

& ,f6.3,3x,f7.4) 

395 continue 

400 continue 

410 goto 5 

end 

C THIS FUNCTION COMPOTES THE INCOMPLETE GAMMA FUNCTION, 
function gamt(a,x) 
if (x.eq.O)goto 40 
if(x.gt. 100) goto 50 
sum=1 ,/a 
an=1.0 
old=sum 

33 old=old*x/(a+an) 

if (old/suffl-1 .6-6)20,10,10 

10 an=an+1 . 

sum= sum+old 
if(an-300.)33,33,12 

12 continue 

20 xxx=(a*alog(x)+alog(sum)-x) 

if (xxx.lt. -80. )goto 40 
gamt=(exp(xxx)) 
goto 60 

40 gamt=0.0 

goto 60 

50 gamt=ganima(a) 

60 return 

end 

c This function computes the gamma function by a Stirling approx. 
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function gamma(y) 

x=y+1 . 

pi=3.14l59 

stirl =!./(! 2. *x) 

stir2=l./(288.*x'^*2.) 

Stirs =~ 1 3 9 . / ( 5 1 840 . *x«»3 . ) 
stir4=-5T1 ./ (2488320. »x»*4.) 
stir=1 +stir1 +stir2+stir3+stir4 
gamma= exp ( -x) *x** ( x- , 5 ) *sqrt ( 2 . *pi ) * st i r/y 
end 

cS************#*****^******.*****^****?*****************#****?:*##*??****** 

subroutine WATCN ( ta , sut , nu , gamsw) 
real nu,nut 

dimension sutt(1 1) ,nut(l 1) ,gamst(l 1) 

data sutt/75. 6, 74. 9, 74. 2, 73.5, 72. 0,72. 1,71.4, 70. 7, 70.0, 

& 69 . 3 , 68 . 6 / 

data nut/17 .93e-3 , 1 5 . 1 8e-3 , 1 3 . 09e-3 , 1 1 .44e-3 , 1 0 . 08e-3 , 8 . 94e-3 , 
& 8.e-3,7.2e-3,6.53e-3,5.97e-3,5.94e~3/ 

data gamst/0. 99987, 0.9 999999, 0.9 9 973, 0.99913, 0.9 9823, 0.997 08 
& ,0.99568,0.99406,0.99225,0.99025,0.98807/ 

if(ta.gt.50.)go to 10 
ita=if ix( ta*. 2)+1 
frac=ta~float(ifix(ta)) 
ita1=ita+l 

sut= ( sut t( ital ) -sut t ( ita) ) *0 . 2*frac+sutt (ita) 

nu= (nut ( ital ) -nut (ita) ) *0 .2*frac+nut ( ita) 

gamsw=( (gamst(ital ) -gamst (ita) ) *. 2*frac+gamst(ita) ) *980 . 

return 

10 sut=sutt(11) 

nu=nut(l 1) 
gamsw= gainst (11) 
return 
end 
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c DIMENSIONLESS EXFILTRATION DIFFUSIVITY 
function fie(d) 
dimension y(6) 

data y/0. 18, 0.1 1,0. 077,0.056,0.044, 0.034/ 
if(d.gt.7.) goto 10 
x=d-1 , 

if(x.lt.1 .)x=1. 

i=ifix(x) 

frac=x-float(i) 

y1=alog(y(i)) 

y2=alog(y(i+D) 

fie=exp( ( y2~y1 ) *frac+y1 ) 

return 

10 fie=.034 

return 
end 

csssssssss^efsssasswftasss^s^s^fsx^aj^^asssfsees****************************** 
subroutine PR0BZ(mnu,p1 ,prob,ak) 
c THIS PROGRAM COMPUTES THE CDF OF NORMALIZED PRECIPITATION. 
real*8 fac(500) 

real®8 x,a,dlog,gama,gamlid,eps 
realms m,k,w,t,2,zl,zu,inz 
real *8 p1-,mnu 

real*8 xold,xsum,sum1 , sum2, tot, vtot, void, vnew 
integer v,vm,vmax 

if(pistol.eq,2) goto 301 
do 300 j=1,500 
vtot=0,0d0 
do 700 iv=1,j 

700 vtot=vtot+dlog(dble(float(iv)) ) 

fac( j) =vtot 

300 continue 

301 continue 
eps=l .e-5 
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pistol=2 

w=mnu 

t=l. 

k=ak 

c INITIALIZING VALUES 
c 

m=w*fc 

z=p1 

vm=if ix( sngl ( m) ) 
vinax=ifix(sngl(3 .%)) 

3 x=m»k»z 

ii=0 
iO=1 

sum1=0,0d0 
sum2=0.0d0 
13 v=vm-ii 

if(v.eq.O)goto 500 
23 if (v.eq.vmaxjgoto 600 

c 
c 

c COMPUTE LOG INCOMPLETE GAMMA DISTRUBITION 

a=dble(float(v))*k 
xold=1.0d0/a 
xsum=1 ,0d0/a 
i=1 

100 xold=(xold/(a+i))*x 

xsum^xsuffl+xold 

if((xold/xsum) .le.eps)goto 200 
i=i+1 
goto 100 
200 continue 

call mlgama(a,gamm,ier) 
gainlid=a*dlog( x) -x+dlog(xsum) -dble( gamm) 
e 

c COMPUTE THE SUMMATION OF ALL V TERMS 
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vold=dble(float(v) ) *dlog(m) -faG(v)+gamlid-in 

if ( void . le . -85 . ) vold=-85 . 

vnew= dexp ( void ) 

if (v.gt.vm)goto 800 

sum1 =sum1 +vnew 

if((vnew/sum1) .le.eps)goto 500 

ii=ii+1 

goto 13 

500 v=vm+jj 

goto 23 

800 sum2=sum2+vnew 

if(( vnew/sum2) .le.eps)goto 600 

goto 500 
c 

e COMPUTE CDF OF NORMALIZED PRECIPITATION 
c 

600 if(m.gt.85.)Di=85. 

prob=sum1 +suia2+dexp(-m) 

return 

end 
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